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This  report  examines  circular  error  probable  (CEP) 
approximation  methodologies.  The  approximations  and 
methodologies  explored  in  this  report  use  statistical 
distributions  and  mathematical  methods  that  may  be 
unfamiliar  to  some  people.  It  is  assumed  that  the  reader 
has  a  general  understanding  of  statistical  distributions, 
and  their  associated  parameters  as  well  as  how  to  manipulate 
these  distributions. 

This  thesis  was  suggested  by  persons  from  the  4220th 
Weapon  System  Evaluation  Squadron,  Strategic  Air  Command  and 
should  be  useful  to  all  offices  interested  in  calculating 
CEPs  for  missiles,  rockets,  bombs,  bullets,  etc. 
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Several  approximation  techniques  are  currently  used  to 
estimate  CEP.  These  techniques  are  statistical  in  nature, 
based  on  the  bivariate  normal  distribution  of  the  crossrange 
and  downrange  miss  distances  of  sample  impacts  of  weapon 
systems.  This  thesis  examines  four  of  the  most  widely  used 
approximation  techniques  (Grubbs-Pat  n aik/chi-square,  Grubbs- 
Pa tna ik/ W i 1 s on-H 1 1 f e r ty ,  modified  RAND-234,  and  correlated 
bivariate  normal),  compares  their  results  with  the  results 
and  computational  effort  required  by  established  numerical 
integration  techniques,  determines  the  relative  accuracy  of 
each  technique  in  various  regimes  of  the  b i as / e 1 1  i  p t  i  c  i  t y 
parameter  space.  Included  in  this  report  is  a  tutorial  on 
the  subject  of  CEP  meant  to  serve  as  a  general  introduction 
to  how  to  calculate  CEPs  with  some  of  the  popular 
approximation  techniques. 

In  general  it  was  found  that  each  of  the  approximation 
techniques  is  best  in  some  regime  of  the  parameter  space 
with  the  G r ubb s -Pa t na i k / c h i - s q ua r e  technique  being  the  most 
reliable  estimator.  For  fast  calculations  of  CEP,  the 
correlated  bivariate  normal  and  the  "exact”  method  may  not 
be  feasible  because  both  are  computationally  rigorous  and 
require  from  2  minutes  to  several  hours  of  computer  time  (on 
a  personal  computer)  to  give  an  estimate  of  CEP. 


EXAMINATION  OF  CIRCULAR  ERROR  PROBABLE 
APPROXIMATION  TECHNIQUES 

I  .  Introduction 

Background 

Circular  error  probable  (CEP)  Is  defined  as  "...  the 
radius  of  a  circle  centered  at  the  target  or  mean  point  of 
impact  within  which  the  probability  of  impact  is  0.5" 
(16:1).  That  is: 


where 


// 

CCEP 


F  (  x  ,  y  )  dx  dy  =  0.5 


(1.1) 


C C e P  =  a  circle  centered  at  the  target  with  radius  CEP 
F ( x , y )  is  the  bivariate  normal  density  function 
x  Is  the  crossrange  miss  distance 
y  is  the  downrange  miss  distance. 


Given  this  definition,  circular  error  probable  (CEP)  is  used 
to  measure  a  weapon  system's  impact  accuracy.  The  CEP  is 
determined  by  taking  test  impact  data  and  approximating  the 
actual  CEP  of  the  missile  system.  Of  course,  the  more  data 
points  -  or  weapon  system  tests  -  collected,  the  better  the 
approximation  of  CEP. 

There  are  several  ways  to  approximate  CEP.  One  group 


of  techniques  is  non -p a r a m e t r  i  c  approximation  methods. 

These  approximation  techniques  make  no  assumptions  about  the 


underlying  distributions  of  the  impacts.  One  non-parametric 
method  simply  uses  the  sample  median  as  an  estimate  of  CEP. 
Non-parametric  estimations  are  only  useful  in  cases  where 
one  has  available  a  large  number  (  greater  than  30  )  of 
sample  impacts.  When  flight  testing  missiles  it  is  not 
practical  to  expect  to  have  a  large  number  of  data  points 
(usually  no  more  than  15  are  available)  (2),  so  non- 
parametric  tests  are  of  limited  use  and  are  not  presented  in 
this  study. 

The  most  common  parametric  methods  for  estimating  CEP 
are  based  on  the  assumption  that  the  impacts  are  normally 
distributed:  closed  form  integration  of  the  bivariate 

normal  function,  numerical  integration  of  the  bivariate 
normal,  or  algebraic  approximations  of  the  bivariate  normal 
function  (16:2-4).  Closed  form  integration  of  the  bivariate 
normal  distribution  is  of  limited  use  because  it  can  only  be 
accomplished  for  the  case  of  n o n -c o r r e 1  a t e d  samples  with 
means  of  zero  and  equal  standard  deviations  (16:2-3). 

On  the  other  hand,  numerical  integration  techniques  are 
useful  in  estimating  correlated  samples  with  non-zero  means 
and  unequal  variances.  Numerical  integration  techniques 
yield  a  probability  of  impact  given  a  circle  of  known 
radius.  However  in  approximating  CEP  it  is  of  interest  to 
derive  a  radius  given  a  known  probability  (i.e.  0.5). 

Therefore  >ne  must  come  up  with  an  initial  estimate  of  CEP 
and  then  iterate  the  numerical  integration  until  a 


probability  of  0.5  is  achieved  (3). 


It  is  important  to  note  that  numerical  integration 
techniques  for  estimating  CEP  are  in  fact  estimates  of  CEP. 
They  are  the  most  accurate  estimates  available  under  the 
assumption  of  normality.  Because  of  the  fact  that  they  are 
most  accurate  they  are  often  referred  to  as  "exact" 
solutions.  Although  the  use  of  an  "exact"  method  provides 
the  most  accurate  results,  there  are  advantages  and 
disadvantages  to  its  use. 


While  numerical  integration  methods  provide  good 
estimates  of  CEP,  and  can  be  used  to  evaluate  the 
accuracy  of  the  other  CEP  approximation  methods,  they 
require  considerable  computer  time,  and  are  usually 
impractical  in  flight  test  analysis  (16:3). 


The  key  point  is  that  numerical  integration  can  be  used  as  a 
standard  to  evaluate  the  accuracy  of  other  CEP  approximation 
methods. 

The  preferred  method  of  approximating  CEP  using 
numerical  integration  of  the  bivariate  normal  density 
function  is  infinite  series  expansion.  When  calculating 
CEP,  one  is  solving  for  the  radius  within  which  the 
probability  of  hit  is  0.50. 


POO 

where, 


d  xd  y  (1.2) 


•//— 

JJ  2*SxSy 


(x,  y)  =  downrange  and  crossrange  miss  distances; 
(x,  y)  =  downrange  and  crossrange  sample  means; 


(SX)  Sv)  =  downrange  and  crossrange  sample  standard 
deviations  . 

When  P(R)  is  equal  to  0.50,  then  R  is  the  estimate  of  the 
CEP.  Smith  suggests  an  iterative  approach  to  solve  for  CEP 
(17:1-6).  One  approximates  this  integral  (1.2)  with  an 
infinite  series  and  expands  it  until  the  inner  terms 
approach  zero.  This  method  of  approximating  CEP  is  referred 
to  as  the  "exact"  method.  The  calculations  required  by  the 
"exact"  method  can  be  extremely  rigorous  as  will  be 
illustrated  in  Chapter  III.  The  number  of  iterations  for 
the  series  to  converge  can  use  much  computer  time. 

Another  numerical  integration  technique  in  current  use 
is  the  correlated  bivariate  normal  method  (CBN).  The  CBN  is 

sometimes  a  faster  approximation  than  the  "exact"  method, 
however,  it  also  requires  considerable  time  to  converge  to 
the  CEP. 

Since  numerical  integration  techniques  may  be 
impractical,  it  is  essential  to  have  fast  approximation 
techniques  that  are  fully  tested  and  validated.  The 
algebraic  approximations  give  fast  results  and  are  accurate 
over  certain  regimes  of  the  parameter  space  (Three 
parameters  can  be  used  to  characterize  the  probability 
distributions  of  impacts:  bias,  ellipticity,  and 

correlation.  It  will  be  shown  that  correlation  can  be 
removed  from  sample  data  so  that  one  is  left  with  two 
parameters).  Many  different  algebraic  approximations  of  CEP 


exist.  Three  of  these  algebraic  approximation  methods  are 

of  particular  interest  to  this  study  of  CEP  methodologies: 

Grubbs-Patnaik/chi-squared 
Modified  RAND-234 
Grubbs-Patnaik/Wilson-Hilferty 

These  three  methods  were  suggested  by  the  sponsors  of  this 

study  (4220th  Weapon  System  Evaluation  Squadron,  Strategic 

Air  Command)  as  the  most  frequent  methods  used.  No  work  has 

been  accomplished  to  fully  compare  these  approximations  with 

each  other  and  with  the  numerical  integration  methods. 

Problem  Statement 

Given  non-c o r r e la t e d  sample  impacts,  how  do  common  CEP 
approximation  techniques  (Grubbs-Patnaik/chi-squared, 
Modified  RAND-234,  Grubbs-Patnaik/Wilson-Hilferty,  or  CEN) 
compare  in  accuracy  and  computational  effort  (measured  by 
computer  time)  to  the  "exact”  method  (numerical  integration) 
over  the  possible  range  of  the  parameters  bias  and 
ellipticity. 

Purpose  o f  the  Study 

It  has  been  suggested  that  the  "exact"  numerical 
integration  method  of  approximating  circular  error  probable 
(CEP)  is  not  practical  for  flight  test  analysis  because  of 
the  inordinate  amount  of  computer  time  it  takes  for  this 
method  to  produce  an  approximation  of  CEP  (16:3).  This 
study  examines  four  methods  of  approximating  CEP  and 
compares  them  to  the  "exact"  method. 

While  examining  these  approximation  methods  this  study 


1-5 


determines  which  method  gives  the  "best"  estimate  of  CEP 
when  compared  to  the  "exact"  methodin  various  regimes  of  the 
b i as / e 1 1 ip t i c i ty  parameter  space.  Additionally,  a  tutorial 
on  CEP  and  the  various  approximation  techniques  for 
calculating  CEP  is  presented. 

The  general  approach  to  this  thesis  is  to: 

--  Apply  the  three  algebraic  approximation  techniques 
and  the  CBN  over  a  wide  range  of  the  bias/ 
ellipticity  parameter  space; 

--  Compare  accuracies  and  computational  effort 

(computer  time)  of  each  of  the  four  approximation 
techniques  with  the  "exact"  numerical  integration 
method; 

--  Analyze  results  to  determine: 

-  which  technique  is  most  accurate  (compared  to 
the  "exact"  method)  in  a  given  regime  of  the 

b ias / e 1 1 ip t ic i ty  parameter  space; 

-  where  any  of  the  techniques  may  fail  to  give  an 
accurate  estimate  of  CEP. 

--  Use  regression  analysis  to  estimate  the  error 

generated  by  each  approximation  technique  over  the 
b ia s / e 1 1 i p t i c i t y  parameter  space  and  add  a 
correction  factor  to  the  calculations. 

This  study  does  not  use  actual  test  data  for  validating 
the  approximation  techniques.  Non-dimensional  values  of 
bias  and  ellipticity  are  used.  Therefore,  correlation  is 


assumed  to  be  equal  to  zero.  Given  actual  test  data,  the 
correlation  can  be  removed  from  the  data  by  translating  the 
data  to  principal  axes  (See  Appendix  B  for  details  on  how  to 
remove  correlation  from  sample  data)  (10:20-21).  The  units 
(feet,  meters,  etc.)  for  CEP  or  any  of  the  parameters  do  not 
affect  the  calculation  of  CEP,  as  long  as  one  is  consistent 
with  the  units  used. 

Seq  uence  o  f  Presentation 

Chapter  II  presents  a  review  of  related  literature  in 
this  field.  Chapter  III  describes  the  three  algebraic 
approximations,  the  CBN,  and  the  "exact"  method.  Also 
Chapter  III  explains  the  data  collection  process,  computer 

program  and  the  regression  analysis  used  to  estimate  the 
errors  between  the  algebraic  approximations  and  the  "exact" 
solution. 

Chapter  IV  discusses  the  results  of  this  study  in  terms 
of  which  approximations  of  CEP  are  "best,"  and  in  what 
regimes  of  the  b ia s / e 1 1 i p t i c i t y  parameter  space  each  should 
be  used.  The  results  of  the  regression  analysis  are  also 
discussed  in  Chapter  IV.  Chapter  V  summarizes  the  findings 
of  this  study  and  makes  recommendations  for  the 
implementation  of  these  results  and  for  further  research. 

The  appendices  present  the  tutorial  on  circular  error 
probable,  the  derivation  of  some  of  the  mathematical 
formulae  used  in  this  study,  computer  program  listings,  and 
some  sample  output.  Appendix  A  contains  the  tutorial  on 


II.  Literature  Review 


Past  work  in  the  area  of  approximating  circular  error 
probable  can  be  divided  into  two  groups:  theoretical  work 
and  applied  work. 


Theore  tical  Work 


The  theoretical  work  in  the  area  of  estimating  CEP  is 

concerned  primarily  with  deriving  the  various  approximation 

techniques.  Chapter  III  will  give  the  details  about  the 

derivations  of  each  approximation  method. 

The  theoretical  groundwork  for  two  of  the  approximation 

techniques  of  interest  to  this  study  was  set  by  Frank  E. 

Grubbs.  He  combined  the  work  of  Patnaik,  Wilson,  and 

Hilferty  to  formulate  two  methods  of  approximating  CEP. 

These  are  referred  to  as  the  G r ubb s -Pa tna i k / c h i -s q ua r e ,  and 

the  G rub bs -Pa tna ik/ W i Is  on  Hilferty  (7).  Grubbs  defined  the 

problem  of  interest  as: 

that  of  finding  the  probability  of  hitting  a 
circular  target  ...  whether  the  delivery  errors 
are  equal  or  unequal  and  also  for  point  of  aim  or 
center  of  impact  of  the  rounds  either  coinciding 
with  the  target  centroid  or  offset  from  it.  ...  It 
therefore  appears  desirable  to  record  a  straight¬ 
forward,  unique,  and  rather  simple  technique  for 
approximating  probabilities  of  hitting  for  all  of 
the  various  cases  referred  to  above  (7:51). 

Grubbs'  methods  are  being  used  by  several  organizations  to 

calculate  CEP:  Strategic  Air  Command  (missiles  and  bombers) 

(4),  the  Army  Missile  Command  (18),  and  USAF  Foreign 

Technology  Division  to  name  a  few. 


In  addition  to  Grubbs'  two  approximation  methods  the 


Rand-234  method  is  presented.  This  method  is  .iot  based  on 

Grubbs'  work.  The  RAND-234  method  of  approximating  CE?  is 

simply  a  mathematical  expression  derived  to  approximate 

RAND-234  tables  of  probabilities  (14).  The  RAND-234  tables 

were  contained  in  RAND  Report  R-234. 

R-234  tables  contain  the  probabilities  of  missing  a 
circular  target  of  a  given  radius  with  a  weapon  system 
of  known  accuracy  aimed  at  some  point  offset  from  the 
center  of  the  target.  R-234  assumes  that  the  weapon 
system  accuracy  can  be  described  by  a  radially 
symmetric  Gaussian  distribution  (14:2). 

A  least  squares  regression  was  used  to  derive  this  formula. 

The  reason  that  this  method  was  developed  was  that  "...  a 

formula  facilitates  computer  calculation  of  CEP  and  obviates 

the  need  to  look  up  CEPs  in  the  R-234  tables"  (14:2).  This 

regression  resulted  in  a  third  order  polynomial  for 

calculating  CEP. 

Two  other  CEP  approximation  techniques  developed  in  the 
literature  are  numerical  integration  techniques.  One  is 
called  the  correlated  bivariate  normal  (CBN),  and  the  other 
is  referred  to  as  the  "exact"  method  because  it  is 
considered  the  most  accurate  approximation  (it  is  also  the 
most  computationally  rigorous  of  the  approximations). 

The  CBN  was  developed  "...  to  calculate  CEP  about  the 
target  point  given  that  the  center  of  the  distribution  of 
impacts  is  located  at  some  distance  (  r a n g e - t r a ck )  defined  by 
an  impact  bias  vector"  (15:1).  The  CBN  solves  for  CEP  by 
expressing  the  bivariate  normal  distribution  in  polar 
coordinates  and  then  integrating  over  r  and  summing  over  9. 


The  resulting  expression  must  then  be  solved  iteratively  for 
the  CEP  value  (r)  that  results  in  a  probability  of  0.5 
(15:1-5).  Typically  40  steps  were  used  in  this  integration. 
Jones,  however,  notes  the  accuracy  of  the  CBN  can  be 
improved  by  increasing  the  number  of  steps  used  in  the 
integration  (11:1). 

Another  numerical  integration  technique  for 
approximating  CEP  is  the  "exact"  method.  The  "exact"  method 
is  described  in  detail  by  Smith.  This  method  involves 
Taylor  series  and  infinite  series  expansions  to  estimate  a 
bivariate  normal  distribution  with  zero  correlation.  After 
expansion,  one  must  solve  for  the  radius  (r)  within  which 
the  probability  of  hit  is  0.5.  This  is  done  by  iterating  on 
r  until  the  series  expansion  converges  to  0.5.  Whereas  the 
CBN  only  involved  summation  over  one  variable  (9),  the 
"exact"  method  sums  over  both  r  and  9  (16). 

Applied  Work 

Comparisons  of  some  of  these  methods  have  been 
accomplished,  but  never  has  there  been  a  study  to  compare 
all  five  of  these  approximation  methods  at  the  same  time. 

Smith  compared  the  RAND-234,  G r u bb s -Pa t na i k / c h i -s q ua r e , 
and  the  G r u bb s -Pa t n a i k / W i 1 s on -H i 1 f e r t y  methods  with  the 
Infinite  series,  "exact"  approximation  (16).  In  this  work 
Smith  showed  that  one  can  produce  a  very  significant  error 
in  calculating  CEP  by  neglecting  the  rotation  to  principal 
axes  in  the  case  of  correlated  samples.  This  error  was  as 


much  as  2  3  Z  In  one  case.  Also  Smith  varied  eccentricity 


while  holding  bias  constant  and  vice  versa.  She  did  not, 
however  show  the  effect  of  varying  both  bias  and 
eccentrlcity  at  the  same  time  (15:9). 

Other  works  have  also  compared  approximation  methods. 
Jones  compared  the  RAND-234,  CBN  (with  40  and  400 
iterations),  Grubbs-Patnaik/chi-square,  and  the  "exact." 
Table  2-1  shows  the  values  of  ellipticity  and  bias  used  in 
the  16  cases  examined  by  Jones. 

Table  2-1 
Cases  (11:3) 


(Note:  the  units  in  this  analysis  do  not  matter  as  long  as 

one  is  consistent  with  the  units  --  once  in  meters,  always 
in  meters,  etc.) 

Table  2-2  shows  the  CEPs  calculated  for  each  of  the 
techniques  (RAND-234,  CBN  with  40  steps,  CBN  with  400  steps 


Grubbs-Patnaik/chi-square,  and  "exact")  for  each  of  the  16 
cases  In  table  2.1.  Table  2.2  Illustrates  that  the  "best" 
approximation  varies  depending  on  the  values  of  ellipticity 
and  bias.  However,  this  paper  did  not  explicitly  delineate 
where  one  method  becomes  more  inaccurate  than  another. 

Tahle  2-2 
Results  (11:4) 


Case 

RAND 

CBN( 40) 

CBN(400) 

Grubbs 

Exact 

1 

801.18 

805.98 

818.81 

824.41 

820 .27 

2 

830.73 

833.54 

845.20 

847.12 

846.56 

3 

865.21 

866.81 

877.67 

876.89 

878.97 

4 

904 .28 

905.23 

915.57 

913.48 

916.82 

5 

947.60 

947 . 08 

957.09 

954.52 

958.34 

6 

994.80 

992  .  5  1 

1002 .37 

999.71 

1003  .6  1 

7 

1045.55 

1041 .44 

1051  .28 

1048.37 

1052.55 

8 

1099.48 

1093.69 

1  103.67 

1100.02 

1104.97 

9 

1156.25 

1149.01 

1159.28 

1154.41 

1160.62 

10 

1215.51 

1207.01 

1217.72 

1211.11 

1219.12 

1 1 

1276.90 

1267 . 24 

1278.54 

1269.88 

1280.01 

12 

1340 . 08 

1329.22 

1341.24 

1330.47 

1342.80 

13 

1404 . 69 

1392.53 

1405.40 

1392.65 

1407.04 

14 

1470 .39 

1458 .07 

1471 . 90 

1473.64 

1457.47 

15 

1536.82 

1521.86 

1536.68 

1521.02 

1538.52 

16 

1603.64 

1587.43 

1603.32 

1586.85 

1605.27 

By  calculating  the  error  of  each  approximation  relative  to 
the  "exact"  method,  one  can  determine  which  approximation  is 
"best"  for  each  case.  Equation  2.1  is  used  to  determine  the 
relative  error  (RE)  (16:9). 


RE  = 


CEPapprox  ~  CEPexact 


(2.1) 


CEP 


exact 


For  example,  table  2-3  contains  relative  errors  for  five  of 
the  example  cases. 


2-5 


Relative  Errors 


Case 

RAND 

CBN{ 40 ) 

CBN (400) 

Grubbs 

1 

-0.023 

-0.017 

-0.002 

0.005 

5 

-0.011 

-0.012 

-0.001 

-0.004 

10 

-0.003 

-0.010 

-0.001 

-0.007 

13 

-0.002 

-0.010 

-0.001 

-0.010 

16 

-0.001 

-0.011 

-0.001 

-0.012 

For  all  cases  the  CBN  with  400  steps  is  the  "best."  However 
as  stated  before  this  method  is  computationally  rigorous 
much  like  the  "exact"  method.  So  if  one  needs  a  "quick" 
estimate  for  case  1  the  Grubbs-Patnaik/chi-square  is  the 
"best."  However,  by  case  10  the  RAND-234  gives  the  "best” 
estimate  of  the  fast  methods,  and  by  case  16  the  RAND-234  is 
as  accurate  as  the  CBN  with  400  steps. 

As  mentioned  before,  no  study  has  been  accomplished  to 
analyze  the  accuracy  of  each  of  the  five  methods  here 
presented  over  as  wide  a  range  of  values  as  those  in  this 


studv. 


III.  Methodology 


This  chapter  presents  the  methodology  used  to  analyze 


the  CEP  approximation  techniques.  The  variables  used  in  the 


methods  (mean,  standard  deviation,  and  correlation)  are 


introduced,  followed  by  a  presentation  of  the  mathematics  of 


each  of  the  approximation  methods.  Next  the  framework  of 


the  attempted  regression  analysis  is  discussed,  followed  by 


a  description  of  the  process  by  which  data  was  generated  for 


this  study.  Finally,  the  interactive  computer  program 


included  in  appendix  E  is  described. 


The  Variables 


Circular  error  probable  can  be  expressed  as  a  function 


of  the  standard  deviation,  mean,  and  correlation  of  the 


downrange  and  crossrange  miss  distances  of  a  set  of  sample 


impacts.  Correlation  can  be  removed  from  sample  data  by 


performing  a  simple  rotation  of  axes  into  a  non-co r re  la  ted 


coordinate  system  (see  Appendix  B  for  the  this  rotation 


technique  and  the  formula  for  calculating  correlation) 


(10:21).  This  study  assumes  one  is  using  the  standard 


deviations  and  means  from  uncorrelated  samples. 


The  unbiased  estimators  of  standard  deviation  and  mean 


are  used  here.  Equations  3.1  and  3.2  show  the  formulae  used 


in  calculating  the  down  range  and  cross  range  sample  means: 


x  -  l  x  i  / n 

i  =  1 


(3.1) 


*  A*  S.‘\*  '*• 


where  x  refers  to  cross  range  miss  distance,  y  refers  to 


down  range  miss  distance,  and  n  is  the  total  number  of 
sample  impacts.  Equations  3.3  and  3.4  are  the  formulae  for 
calculating  the  sample  standard  deviations: 


S 


x 


/ 


n 


n 


I 

i  =  1 


(yi2 


nx2  )  /  (  n 


ny  2  )  /  (  n 


1) 


1  ) 


(3.3) 

(3.4) 


From  these  variables  the  parameters  of  bias  and  ellipticity 
are  obtained. 

bias  =  V  *2  +  y2  (3.5) 

ellipticity  =  Sx/Sy  (3.6) 

The  relative  errors  between  the  approximation  techniques  and 
the  "exact"  approximation  vary  depending  on  the  bias  and 
ellipticity  of  a  given  set  of  sample  impacts. 


Grubbs-Patnaik/Chl-Square 

The  G r u bb s -Pa tn a i k / c h i -s q ua r e  approximation  method  is 
based  on  the  work  of  Frank  E.  Grubbs  (7).  This  method  uses 
the  fact  that  the  bias  is  a  sum  of  noncentral  chi-squares 
and  solves  for  the  radius  that  gives  a  probability  of  0.5 
(i.e.  P[(x-x)2/Sx2  +  (y-y)2/Sy2  <=  R2]  =  0.5)  (7:55). 


This  is  the  approximation: 


CEP  =  V  kv/2m  (3.7) 

where 


m  = 

(Sx2  +  s  y  2  +  + 

y2> 

(3 

.8) 

v  = 

2(SX4  +  2p2Sx2Sy2 

+  s  y  4 )  + 

4(x2Sx2  + 

2xypSxSy  +  y  2  S  y  2  ) 

(  3 

.9  ) 

k  = 

F“ 1 ( 0 . 5  ) 

(  3  . 

10) 

df  • 

=  2  m  2  /  v 

(3. 

ID 

and  F  is  the  chi-square  distribution  function  with  df 
degrees  of  freedom  (1:5).  In  this  formula  p  is  the 
correlation  coefficient.  If  the  rotation  described  in 
appendix  B  is  performed  then  the  terms  with  p  in  them  can  be 
Ignored  (p  =  0). 

In  the  computer  programs  used  In  this  study  the  inverse 
chi-square  function  was  estimated  using  table  values  found 
in  the  CRC  Standard  Math  Tables  (5:547).  This  table  gives 
values  of  F“*  for  integer  values  of  df,  but  F“*  is  a 
continuous  function.  To  obtain  values  of  F-^  for  non-Intege 
values  of  df  a  simple  linear  interpolation  was  used  (Note: 
F“1  is  not  a  linear  function,  but  for  the  sake  of  simplicity 
the  linear  approximation  was  used  for  the  interpolation. 
Comparing  the  results  obtained  using  this  simple  linear 
Interpolation  and  results  of  the  G r u b b s - Pa t n a i k / c h i -s q u a r e 


from  other  studies  that  used  a  more  exact  expression  for  F  “  * 
showed  that  the  was  little  or  no  loss  of  accuracy.). 

Grubbs-Pa  tnalk/mison-Hllferty 

The  G r u b b s - Pa t na i k / W i 1 s on - H i 1 f e r t y  approximation 
technique  was  developed  as  a  modification  of  the  Grubbs- 
Patnaik/chi-square  method.  The  Grubbs-Patnaik/Wilson- 
Hilferty  method  transforms  the  chi-square  to  approximate 
normal  variables.  This  method  does  not  use  the  chi-square 
function  described  in  the  previous  section.  The  expression 
for  CEP  used  in  the  Grubbs-Patnaik/Vilson-Hilferty  method 
i  s  : 

CEP  =  ^  m { 1  -  [ v/ ( 9m2  )  ]  3)  (3.12) 

where  m  and  v  are  as  defined  in  the  G r u b b s - P a tna Ik / c h  i  - 
square  method  (equations  3.8  and  3.9)  (1:6). 

Modified  RAND-234 

The  Modified  RAND-234  method  is  a  fit  of  a  cubic 
polynomial  to  a  table  of  CEP  values.  Pesapane  and  Irvine 
used  regression  analysis  to  "...derive  a  mathematical 
expression  for  circular  error  probable  (CEP)  which 
approximates  probability  tables  contained  in  RAND  Report  R- 
234”  (14:2).  This  is  the  Modified  RAND-234  method: 

CEP  =  CEPMpT( ’ .0039  -  0.0528v  +  0.4786v2  -  0.0793v3)  (3.13) 

where  CEPMpj  is  the  CEP  centered  on  the  mean  point  of  impact 


(this  CEP  must  be  translated  to  a  CEP  centered  on  the 


3-5 


normal )  : 


? ( r*  )  = 


2  IT  S  XS  y  V  1"P' 


2  tt  r 

// 

0  0 


exp f -( a r 2+2 b r+c ) ]  rdrd©  (3.19) 


where  r  *  ( P  ( r  * )  =  P(CEP)  =  0.5)  is  the  radius  for  which  one 
is  solving,  a  and  b  are  functions  of  the  change  in  9 
(eauations  3.20  and  3.21),  and  c  is  a  constant  (equation 
3.22).  Finally,  equations  3.23  and  3.24  give  the  formulae 
for  the  CBN  approximation  technique.  In  these  equations  N 
is  the  number  of  intervals  the  integral  is  divided  into  to 
give  this  approximation,  AO  =  2tt/N,  a1  and  bj  are  as  given 
in  equations  3.20  and  3.21  for  the  current  value  of  0 
(  2  tt  i  /  N )  ,  and  c  is  as  given  in  equation  3.22.  The  more 
intervals,  the  more  accurate  the  approximation  (15:1-5). 


a  (  9 )  = 


2  (  1 


1  /sin 

-P2)  \  sx 


9  9  > 

^0  2psinOcosO  cos  0 


(3.20) 


S  XS  y 


b  (  9  )  = 


2  (  1 


x  s  i  n  0  PxcosO  +  P  y  s i n  0  v  c  o  s  0 


SxSy 


) 


(3.21) 


c  = 


2(1 


1  /  x2  2  Px  y  y2\ 
-P2)  ySx2  SxSy  Sy  J 


(3.22) 


erf(z)  =  -ssr  I  exp|-x2]dx 


r 


(3.23) 


To  obtain  a  CEP  using  the  correlated  bivariate  normal 
method  one  must  first  make  a  guess  at  what  the  CEP  actually 
is,  solve  equation  3.23  and  compare  the  value  to  0.5.  If 
the  value  obtained  is  "close  enough"  to  0.5,  one  has  a  CEP. 
However  if  the  probability  given  by  the  CBN  is  not  "close 
enough"  to  0.5  then  one  must  iterate  to  a  value  that  is 
closer  to  0.5.  To  decide  how  close  is  "close  enough"  one 
must  decide  how  accurate  your  estimate  of  CEP  is  to  be.  The 
secant  method  (described  in  Appendix  C)  is  used  to  Iterate 
for  the  value  of  the  CEP  In  the  CBN. 

"Exact"  Method 

Like  the  correlated  bivariate  normal  technique,  the 
"exact"  method  estimates  the  bivariate  normal  function.  In 
order  to  use  this  method  the  two  coordinates  of  miss  must  be 
uncorrelated.  This  method  uses  Taylor  series  and  infinite 
series  expansion  to  estimate  CEP.  The  "exact"  method  takes 
the  bivariate  normal  function  in  polar  form: 


I  f  1  (rco 

P  ( R )  =  J  J  A  a  x  p -  — 


s0  -  x)z  ( r s  in0  -  v ) 


rdrd8  (3.25) 


where 


(3.26) 


2  7T  S  y  S  y 


After  series  expansion  this  is  the  equation  for  the  "exact' 
approximation  technique: 


R  2  oo  co  (k  +  l)!(j  +  l)i 

P  ( R )  =  0.5=  —  D  l  l  xkyi  - 

2  k  =  0  j  =0  (k+j  +  1)! 


(3.27) 


where 


S  xSy 


(2k)! 


exp 


"  l  - - -  (^) 

\8SX2/  1=0  (k-1) ! (21)  !  \S  2  / 


-9?2v  1 


(k+1 ) ! ( k ! ) 2 


(  2  j  )  ! 

j  (  j  + 1 )  •  (  j  !  )  2 


(3.28) 


(3.29) 


ta  ‘  J 


r.2v  i 


=0  (j-i)!(2i) 


•(-) 


(3.30) 


As  with  the  CBN  to  solve  this  one  must  first  make  a 
ouess  at  the  CEP,  and  then  iterate  to  arrive  at  an  answer 
that  is  accurate  enough  for  one's  purpose.  The  secant 
method  was  used  in  this  case  also  (17:1-5). 


simply  divides  the  estimated  CEP  by  one  plus  the  relative 
error.  Recalling  equation  2.1: 


CE^approx  “  ^'’exact 


=  RE 


CEP 


exact 


therefore  , 


(2.1) 


CEP 


exact 


CEP 

approx 
(1  +  RE) 


(3.34) 


The  results  of  this  least  squares  multiple  linear  regression 
are  discussed  in  chapter  IV. 

The  Compu ter  Prog  ram 

The  computer  program  included  in  this  study  is  written 
in  the  Pascal  programming  language  and  it  is  compatible  with 
personal  computers  running  Turbo  Pascal  by  Eorland 
International  (with  some  minor  changes  it  will  run  on  other 
versions  of  Pascal).  The  program  operates  interactively. 

It  will  compute  CEPs  for  each  approximation  technique,  only 
those  requested  by  the  user,  or  It  will  decide  which 
technique  is  the  most  accurate  (other  than  the  "exact" 
method)  for  the  given  parameters  and  return  a  CEP  for  that 
technique  only.  The  program  decides  which  method  is  the 
"test”  using  a  decision  criteria  based  on  the  results  of 
this  study. 

The  computer  runs  to  collect  the  data  for  this  study 


were  performed  on  the  VAX  computer  at  the  Air  Force 


Institute  of  Technology,  and  on  the  SANYO  MBC-550  personal 
computer  (operating  system:  MS-DOS)  owned  by  the  author. 

This  program  is  intended  to  run  on  a  personal  computer  in 
the  offices  of  persons  concerned  with  calculating  CEPs.  The 
programs  were  also  run  on  a  Zenith  Z-150  at  the  Air  Force 
Institute  of  Technology.  There  were  no  differences  in  how 
the  program  ran  on  the  SANYO  and  the  Zenith.  The  timing 
analysis  was  based  on  time  to  run  on  the  SANYO. 

The  program  can  easily  be  adapted  for  use  with  any 
version  of  the  Pascal  language.  The  listing  of  the  computer 
program  is  included  in  appendix  E. 

Data  Collection 

To  collect  data  for  this  study  CEPs  for  all  five  of 
the  approximation  methods  were  calculated  over  a  wide  range 
of  values  of  ellipticity  and  bias.  Since  the  bias  is  the 
radial  distance  of  the  mean  point  of  impact  from  the  target 
centroid,  bias  was  stepped  out  on  the  diagonal  where  x  =  y, 
whereas  ellipticity  was  calculated  by  keeping  Sy  constant 
and  varying  Sx.  The  range  of  e 1 1 i p t i c i t i e s  considered  was 
from  0.05  to  1.0,  and  the  range  of  bias  was  from  50  to  1  5  5  5  . 
Overall,  440  data  points  were  examined  in  this  study.  These 
440  data  points  examined  e  1 1  i  p  1 1  c  I  t  i  e  s  from  1.0  to  0.05  in 
decrements  of  0.05  and  bias  from  1555  to  70  (x  and  y  from 
1100  to  50  in  decrements  of  50),  and  a  CEP  was  calculated 
for  each  of  the  five  approximation  methods  for  all 
combinations  of  ellipticity  and  bias. 


Appendix  D  contains  sample  output  from  the  data 
collected.  The  sample  in  Appendix  D  shows  CEPs  for 
ellipticitv  from  1.0  to  3.1  (in  0.1  decrements)  and  bias 
from  1555  to  141  (x  and  y  from  1100  to  100  in  decrements  of 
100)  . 

Once  all  the  data  runs  were  completed,  the  analysis  of 
the  approximation  methods  could  begin. 


IV.  Analysis  and  Discussion  of  Results 


This  chapter  discusses  the  results  of  the  analysis 
described  in  chapter  III.  Each  of  the  approximation 
techniques  is  presented  with  a  discussion  of  how  accurate 
(compared  to  the  "exact"  method  and  the  other 
approximations)  the  method  is  over  the  range  of  the 
b i a s / e 1 1 i p t  i  c i t y  parameter  space.  Also  considered  in  the 
analysis  of  each  approximation  technique  is  the  computer 
time  each  technique  takes  to  give  an  answer  for  the  CEP. 
Finally  this  chapter  discusses  the  regression  analysis 
developed  to  attempt  to  correct  the  approximations. 

Grubbs-Patnaik/Chl-Square 

The  Grubbs-Patnaik/chi-square  method  gives  relative 
errors  (RE)  that  range  from  -0.0076  to  0.0684  (in  several 
cases  the  absolute  relative  error  was  less  than  0.0001). 
Recall  equation  2.1: 

RE  =  CEPaPProx  ~  CEPexact  (2 

CE^exac  t 

The  Grubbs-Patnaik/chi-square  approximation  technique 
underestimates  the  CEP  28%  of  the  time.  Figure  4-1  shows 
the  portion  of  the  parameter  space  where  the  Grubbs- 
Patnaik/chi-square  method  is  the  best. 

This  approximation  is  the  most  accurate  of  the 
approximation  methods.  It  is  the  most  accurate  method  for 
46%  of  the  data  points.  Additionally  the  Grubbs- 


P  a  t  n  a  i  k /chi-square  method  gives  an  answer  in  an  average  of 
seconds . 
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Grubbs-Patnalk/Wilson-Hilferty 

The  Grubbs-Patnalk/Wilson-Hilferty  method  yields 
relative  errors  ranging  from  -0.054  to  0.0726  (the  smallest 


absolute  relative  error  for  this  method  was  0.0C02).  The 
percentage  of  underestimated  C  E  ?  s  given  by  this  method  is 
10%.  Figure  4-2  shows  where  this  method  is  the  most 
accurate  method.  The  G r u b b s -Pa t na i k / W  i  1 s on -H i 1 f e r t y 
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approximation  is  the  most  accurate  method  for  5%  of  the  data 


points.  Like  the  previous  technique  this  method  gives  an 
answer  in  an  average  of  2  seconds. 

Modified  RAND-234 

Values  of  -0.0125  to  0.0765  (0.0002  is  the  smallest 

Ellipticity 
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absolute  relative  error  for  this  method)  are  the  range  of 
the  relative  errors  given  by  the  modified  RAND-234 
technique.  It  should  also  be  noted  that  the  modified  RAND- 
234  method  underestimates  the  CEP  60%  of  the  time.  Figure 
4-3  shows  where  the  modified  RAND-234  method  is  most 
accurate  in  the  bias/ellipticity  parameter  space.  For  29% 
of  the  data  points,  the  modified  RAND-234  approximation  gave 
the  most  accurate  estimate  of  CEP  when  compared  to  the  the 
other  approximation  methods.  This  approximation  technique 
also  returns  estimates  of  CEP  in  an  average  of  2  seconds. 

Astbury  notes  that  the  modified  RAND-234  method  was 
developed  under  the  boundary  condition  that  ellipticity, 
SS/SL,  is  greater  than  0.25.  However,  the  data  for  this 
study  showed  that  the  modified  RAND-234  is  also  reliable  for 
some  values  of  ellipticity  less  than  0.25. 

Correia  ted  Bivariate  Normal 

The  Correlated  Bivariate  Normal  (CBN)  approximation 
technique  can  be  the  most  accurate  estimate  of  CEP  over  the 
entire  range  of  the  parameter  space  if  one  has  the  time  to 
wait  for  a  result.  This  method  takes  anywhere  from  10 
minutes  to  several  hours  to  produce  a  CEP  that  is  accurate 
to  within  0.001  absolute  relative  error.  The  CBN  can  return 
as  accurate  a  value  as  desired.  The  more  intervals  the 
integral  is  broken  into,  the  more  accurate  the 
approximation.  At  times  the  CBN  takes  longer  than  the 
"exact"  method  to  give  a  CEP.  However,  the  CBN  always 
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converges  to  an  estimate  of  CEP. 

Here  is  an  analysis  of  the  time  it  takes  the  CBM  to 
converge  for  three  different  numbers  cf  intervals  the 
integral  is  broken  into: 


40  intervals 


100  intervals 


10  minutes 


1  hour 


400  intervals  -  16  hours 


As  one  can  see,  the  amount  of  time  it  takes  for  the  CBN  to 
converge  increases  exponentially.  Moreover,  to  achieve 
accuracies  to  within  1%  it  is  often  necessary  to  increase 
the  number  of  intervals  in  the  CBN. 


"Exact" 

The  "exact"  method  is  the  benchmark  against  which  all 
the  other  approximation  techniques  were  measured.  This 
method  takes  anywhere  from  2  minutes  to  2  hours  to  give  an 
answer  for  CEP  when  it  converges.  The  "exact"  method  does 
not  always  converge  to  a  CEP.  For  some  highly  elliptical 
distributions  with  large  biases  where  the  target  centroid  is 
not  within  2  standard  deviations  of  the  mean  point  of 
impact,  the  "exact"  method  diverges.  Figure  4-4  shows  the 
regimes  of  the  b i a s / e  1 1  i  p t  i  c  i  t y  parameter  space  where  the 
"exact"  method  does  not  converge. 

In  the  other  regimes  of  the  parameter  space,  the 
"exact"  will  converge  if  one  has  the  time  necessary  for  it 
to  do  so.  If  quick  results  are  not  necessary,  then  the 
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method  is  the  preferred  method  to  calculate  CEP 


Ellipticity 
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Figure  4-5 


Method  Across  the  B  i  a s  /  E 1 1 i p t i c i t y  Parameter  Space 
/  =  RAND-234 

*  =  G r u bbs -Pa t na ik/c h i -s q ua r e 
-  =  G  r  u  b  b  s  - P  a  t  na  i k  / W  i  1  s  o n-H  i  1  f  e  r  t  y 

=  Area  where  the  "exact”  does  not  converge 


X 


For  some  values  of  ellipticlty  and  bias,  the  "exact"  method 
converges  so  quickly  that  it  may  be  unnecessary  to  use  a 
less  accurate  approximation.  The  folowing  is  a  list  of  how 
long  it  takes  the  "exact”  method  to  give  a  result  for 
various  values  of  the  parameter  space: 

5  -  15  minutes 

for  Sx  >  x  and  Sy  >  y 
20  -  30  minutes 

for  0.8  <  Sx/Sy  <  1.0  and  Sx  <  x  and  Sy  <  y 
30  -  90  minutes 

for  0.5  <  Sx/sy  <  0.8  and  Sx  <  x  and  Sy  <  y 
>  90  minutes 

for  all  other  values  of  ellipticlty  and  bias. 

The  notable  exceptions  to  the  above  analysis  comes  when  Sx 
equals  x  or  Sy  equals  y.  For  these  values  the  "exact" 
approximation  method  gives  an  answer  for  CEP  in  less  than  2 
minutes.  Clearly  there  are  times  when  it  is  feasible  to  use 
the  "exact"  method  for  calculating  CEP.  Again,  use  of  the 
"exact"  method  would  depend  on  the  amount  of  computational 
effort  one  is  willing  to  expend  to  achieve  an  answer. 

Figure  4-5  maps  the  entire  b i a s / e 1 1 i p t i c i t y  parameter 
space  examined  in  this  study  and  where  each  of  the  methods 
give  the  most  accurate  estimates  of  CEP.  It  should  be  noted 
that  in  those  regimes  where  the  "exact"  method  will  not 
converge,  one  can  use  the  CBN  to  achieve  a  highly  accurate 
estimate  of  CEP  (if  the  time  is  available). 


Regression  Analysis 

As  part  of  this  study  a  multiple  linear  regression  was 
performed  on  three  of  the  approximation  techniques  (Grubbs- 
Pat  na  ik  /  c  h  i  -s  q  u  a  re  ,  Gr u b bs -Pa tna i k / W  i  1 s on-H  i  1 f e r t y ,  and  the 
modified  RAND-234).  This  regression  analysis  was  performed 
in  the  hope  of  providing  a  correction  factor  for  the 
approximations  so  they  would  more  closely  estimate  CEP. 

This  regression  analysis  was  based  on  the  relative  errors. 

Here  are  the  three  regression  equations  developed  using 
the  regression  procedure  described  in  Chapter  III  (RE  stands 
for  relative  error): 

Grubbs-Patnaik/chl-square 

RE  =  -0.0002  56  +  0 . 0  0  0  3  9  2  (  E 1 1  i  p  )  +  0 .00  000  1  7  2  ( b  i  a  s  )  (4.1) 

Grubbs-Patnaik/Wllson-Hllferty 

RE  =  -0.000472  +  0 . 00069 3 ( Ell  ip  )  +  0 . 0000039 5 ( b  ias  )  (4.2) 

Modified  RAND-234 

RE  =  -0.000403  +  0 .000498( Ellip)  +  0 . 00 000 3  I  6 ( b  i  a s  )  (4.3) 

Once  an  estimate  of  RE  is  determined,  recall  equation  3.30: 


CEP 


exact 


^^approx 
(1  +  RE) 


(3.30) 


is  used  to  come  up  with  a  better  approximation  to  the 
"exact". 


These  regression  equations  give  good  results  for  CEP 


estimates  that  are  greater  than  the  "exact"  CEP.  In  other 
words,  they  do  improve  the  estimate  of  the  CEP.  However 
they  do  not  capture  the  pattern  that  causes  the 
approximation  methods  to  underestimate  the  CEP.  Using  the 
3 ru bbs -Pa tna  ik/ ch i -sq ua r e  method  here  is  an  example  of  the 
results  of  the  regression: 

ellipticity  =  0.5  bias  =  232.84 

CEP(apt»rox.)  =  907.75  CEP(  "exact")  =  906.50 

RE(actual)  =  0.0014  RE ( r e g r e s s  i  on  )  =  0.0004 

Using  the  regression  technique  described  one  would  end  up 
with  a  CEP  of  907.36,  which  has  a  relative  error  of  0.0009. 
That  is  some  improvement  from  the  original.  However,  when 
the  method  underestimates  the  CEP  this  is  a  typical  result: 

ellipticity  =  0.7  bias  *  848.53 

CEP (approx. )  =  1252.7  CEP(”exact")  =  1254.65 

RE(actual)  =  -0.0016  R E ( reg r e s s i on )  *  0.0015 

After  using  the  described  regression  technique  the  CEP  is 
1250.8.  This  is  a  number  further  from  the  "exact”  CEP  than 
the  original  approximation. 

Similar  results  were  obtained  for  each  of  the 
approximation  techniques  when  a  suitable  correction  factor 
was  desired  and  the  approximation  had  underestimated  the 
CEP.  On  the  average,  the  regression  equation  for  the 
Gr ubb s -Pa tna ik/c h i -s q ua r e  method  improves  the  relative  error 
between  the  approximation  and  the  "exact"  by  0.0005.  The 


regression  for  the  Grubbs-Patnaik/Wilson-Hilferty  results  in 
an  improvement  of  0.0008  on  the  average.  The  modified  RAND- 
234  regression  gives  an  average  improvement  in  relative 
error  of  0.0003.  The  fact  remains  that  the  regressions  do 
not  reliably  predict  the  times  that  the  approximations 
underestimate  the  CEP.  It  seems  as  if  the  relative  errors 
of  the  approximation  techniques  are  not  well  behaved  enough 
to  be  captured  in  a  simple  equation.  Appendix  D  table  D-3 
contains  a  listing  of  the  relative  errors  obtained  for  some 
example  values  of  ellipticity  and  bias. 

The  author  also  tried  several  variations  for  the 
regression  independent  variables.  Squaring  the  ellipticity 
and  bias  was  attempted  as  well  as  taking  the  logarithms  of 
the  variables.  All  attempts  at  providing  a  reliable 
correction  factor  for  the  algebraic  approximation  methods 
were  unsuccessful. 

The  scope  of  this  study  was  by  no  means  completely 
exhaustive.  However,  the  author  is  able  to  make  some 
conclusions  and  recommendations  based  on  the  analysis  of  the 
results  of  this  study. 


V  .  Conclusions  and  Recorcir, sedations 

Conclusions 

Circular  error  probable  (CEP)  is  an  important  measure 
of  a  weapon  system's  accuracy.  As  such  it  is  essential  to 
be  able  to  closely  approximate  this  number.  This  study  has 
examined  five  different  methods  of  estimating  CEP  and  from 
this  examination  the  author  has  reached  several  conclusions 
regarding  CEP  calculations. 

First  and  foremost  among  the  author's  conclusions  is 
that  among  the  approximation  methods  studied,  CEP  is  not 
"well  behaved."  What  is  meant  by  that,  is  that  there  is  no 
easily  recognizable  pattern  among  the  approximations  as  to 
which  is  most  accurate  over  a  given  range  of  the  parameter 
space.  This  is  evidenced  by  the  failure  of  the  regression 
analysis  that  was  attempted  in  this  study.  Had  the 
differences  between  the  algebraic  approximations  and  the 
"exact"  method  been  somewhat  "well  behaved"  the  regression 
would  have  provided  a  much  more  reliable  correction  factor. 
However,  this  was  not  the  case. 

Another  observation  is  that  the  algebraic 
approximations  examined  (G r u bb s -Pa t na i k / c h i -s q ua r e ,  Grubbs- 
Pat  na  ik  /  W  i  1  s  on  -H  i  1  f  e  r  t  y  ,  and  modified  RAND-234)  are 
approximately  99%  accurate  for  most  values  of  ellipticity 
and  bias  when  compared  to  the  "exact”  method.  Notable 
exceptions  to  this  are  where  x  equals  Sx  or  y  equals  Sy. 

The  resulting  values  of  ellipticity  and  bias  cause  the 


accuracies  of  the  algebraic  approximations  to  fall  to 
between  90%  and  98%.  However,  for  these  values  the  "exact” 
method  produces  an  answer  in  seconds.  Highly  elliptical 
distributions  cause  the  greatest  errors  in  all  methods  (if 
ellipticity  is  less  than  0.3,  accuracies  drop  to  between  95% 
and  99%). 

The  correlated  bivariate  normal  (CBN)  method  is  the 
most  consistent  method  of  all.  However,  for  this 
approximation  method  to  deliver  accuracies  better  than  99%, 
the  CBN  often  requires  hours  to  produce  an  estimate  of  the 
CEP.  The  only  time  that  the  CBN  would  be  the  preferred 
method  is  when  extremely  accurate  results  (better  than  99%) 
are  needed  for  highly  elliptical  cases  where  the  "exact” 
method  fails  to  converge.  For  these  results  to  be  accurate 
to  better  than  0.1%,  one  must  increase  the  number  of 
intervals  for  the  integral  approximation  to  over  100,  and 
this  required  several  hours  on  the  author's  SANYO  MBC-550 
personal  computer. 

The  "exact"  method,  in  some  regimes,  is  very  practical 
to  use.  If  Sx  is  greater  than  or  equal  to  x  and  Sy  is 
greater  than  or  equal  to  y  the  "exact”  method  converges  to 
the  CEP  very  fast  (less  than  8  minutes).  The  "exact"  method 
has  one  major  drawback  and  that  is  that  it  requires  a 
computer  to  run  the  iterations  (this  is  also  true  of  the 
CBN).  The  algebraic  approximation  techniques  involve  only  a 
few  equations  that  are  solved  once.  These  could  conceivably 
be  programmed  into  a  hand-held  programmable  calculator.  The 


"exact”  method  requires  series  expansion  until  inner  terms 
approach  zero  and  then  an  iterative  scheme  to  achieve  a 
desired  accuracy  in  the  calculation.  In  the  areas  where  the 
"exact"  method  does  not  converge,  one  must  rely  on  the  other 
approximations  for  estimates  of  CEP. 

Recommendations 

Based  on  the  difficulties  this  author  found  in 
programming  these  approximations,  Pascal  is  probably  not  the 
best  language  for  these  calculations.  Pascal  has  some 
limitations  on  the  size  of  the  real  numbers  it  will  accept. 
Turbo  Pascal  assigns  real  numbers  to  a  6  byte  word  size  with 
no  provision  for  double  precision  reals.  This  limits  the 
range  of  real  numbers  to  IE-38  to  1E+38.  Many  of  the 
calculations  required  in  the  approximation  techniques 
produce  numbers  out  of  the  range  of  those  acceptable  in 
Pascal.  The  author  was  reduced  to  using  logarithms  and 
other  manipulations  to  bypass  these  temporary  large  numbers. 
This  resulted  in  some  inefficiencies  in  the  programs  which 
could  lead  to  longer  processing  time.  The  author  also  found 
this  same  problem  with  the  version  of  Pascal  that  is 
available  on  the  VAX  computer  at  the  Air  Force  Institute  of 
Technology.  The  author  is  aware  of  the  existence  of  other 
versions  of  Pascal  that  can  handle  double  precision  real 
numbers.  This  version  of  Pascal  was  not  available  for  use 
in  this  study. 

There  are  other  computer  languages  that  would  perhaps 


lend  themselves  more  readily  to  the  great  number  of 
calculations  required  In  calculating  CEP.  It  may  be  of  use 
to  those  who  calculate  CEPs  to  have  the  programs  furnished 
with  this  study  translated  into  one  of  these  other 
languages . 

In  addition  to  trying  other  computer  languages  this 
author  is  not  convinced  that  there  is  no  way  to  capture  the 
essence  of  the  errors  in  the  algebraic  approximations  in 
order  to  produce  a  correction  factor  for  these  methods. 

This  portion  of  this  study  warrants  further  efforts. 

The  approximation  techniques  examined  in  this  study 
were  all  based  on  the  assumption  that  the  distribution  of 
impacts  follows  the  bivariate  normal  distribution.  It  would 
be  useful  to  have  robust  techniques  developed  that  could 
provide  good  estimates  of  CEP  (given  small  samples)  that  are 
independent  of  the  distribution  of  the  impacts.  This  would 
obviate  the  need  to  make  any  assumptions  about  the 
distribution  of  the  impacts. 

Finally,  this  author  recommends  that  if  one  needs 
highly  accurate  estimates  of  CEP  and  sufficient  time  is 
available,  the  "exact"  method  should  be  used  for  all  values 
where  the  "exact"  converges  to  a  CEP.  In  those  cases  where 
the  "exact"  does  not  converge,  the  CBN  with  at  least  100 
intervals  should  be  used. 

Otherwise  one  will  be  forced  to  use  one  of  the 
algebraic  approximations.  If  one  is  to  chose  a  single 


approximation  technique  to  use,  the  Grubbs-Patn a  ik/chi- 


square  is  the  most  accurate  for  the  greatest  range  of  values 
of  ellipticity  and  bias.  The  analysis  presented  in  this 
study  will  aid  the  user  in  determining  which  approximation 
is  the  "best"  for  one's  purposes. 


Appendix  A 


Circular  Error  Probable  Tutorial 


Introduction 

Whether  one  is  dealing  with  guided  bombs,  unguided 
bombs,  bullets,  missiles,  rockets,  lasers,  etc.,  for  a 
weapons  systems  planner  to  effectively  plan  the  employment 
of  a  weapon  there  must  be  some  measure  of  that  weapon's 
accuracy.  Circular  error  probable  (CEP)  is  one  measure  of 
accuracy  that  is  frequently  used.  CEP  is  defined  as: 

...  the  radius,  centered  on  the  target,  within 
which  the  probability  of  impact  is  0.5. 

To  determine  a  weapon's  CEP  one  must  test  the  weapon  to 
get  a  collection  of  sample  impacts.  This  collection  of 
sample  impacts  may  be  large  (greater  than  30)  or  it  may  be 
quite  small  (take  the  case  of  a  multi-million  dollar  missile 
-  one  wouldn't  launch  30  or  more  of  them  just  to  test  their 
accuracy).  In  the  case  of  a  large  sample  size,  non- 
parametric  techniques  (i.e.  those  that  make  no  assumptions 
about  the  underlying  statistical  distribution)  are  adequate 
for  estimating  the  CEP.  But  no n -p a r a m e t r  i  c  techniques  are 
subject  to  large  errors  when  dealing  with  small  sample 
sizes.  In  the  case  of  small  sample  sizes,  one  must  resort 
to  parametric  methods  for  approximating  CEP.  These 
parametric  techniques  assume  a  particular  statistical 
distribution  and  estimate  the  various  parameters  of  that 


distribution  (i.e.  mean,  standard  deviation,  correlation, 

etc.). 

This  tutorial  gives  a  brief  introduction  to  these 
various  techniques  for  estimating  CEP  starting  with  non- 
parametric  techniques  followed  by  a  discussion  of  parametric 
techniques.  Finally,  results  of  CEP  calculations  for  10 
sample  impacts  of  a  hypothetical  ICBM  warhead  using  each  of 
the  approximation  techniques  is  presented. 

Won-Parametric  Approximation 

A  good  non-parame tr  ic  approximation  for  CEP  is  the 
sample  median.  To  determine  the  sample  median,  one  must 

rank  order  the  sample  miss  distances  (that  is  the  straight 
line  distance  from  the  actual  impact  point  to  the  target). 

The  sample  median  (X)  is  the  middle  statistic  (or  if  there 
are  an  even  number  of  sample  points,  the  median  is  the 
average  of  the  two  middle  order  statistics  of  the  sample 
points). 

X  =  X((n  +  l)/2)  lf  n  is  odd»  ( A .  1  ) 

X  =  f  x ( n / 2 )  +  x( (n  +  1  )/2  )  ^ /2  if  n  is  even*  (A. 2) 

The  sample  median  is  a  good  statistic  to  use  (if  one  has  a 
large  sample)  because  it  requires  no  assumptions  about  the 
underlying  distribution  of  the  impacts. 

Parametric  Approximations 

Parametric  approximations  assume  that  there  is  an 
underlying  distribution  to  the  impacts.  This  report 


and  the  probability  distribution  function  of  the  downrange 
miss  distances  is: 


exp[-0.5(y/Sv)2] 

f  (  y  )  =  — -  -w  <  y  <  eo  (A. 8) 

>Fsy 

If  it  is  assumed  that  S  =  Sx  =  Sy,  and  that  the  crossrange 
and  downrange  miss  distances  are  independent,  then  the  joint 
distribution  of  the  miss  distances  is: 


f  (x,y) 


exp{-0.5((x2  +  y2)/S2]} 
2  m  S2 


(A. 9) 


Changing  to  polar  coordinates,  the  joint  distribution  is: 


f (r,0) 


r  exp [ -0 . 5( r/S  )  2  ] 
2m  S2 


0  <  r  <  °° 

(A. 10) 

0  <  0  <  2n 


To  solve  for  a  particular  radius  r*  such  that  The 
probability  that  r  is  less  than  r*  equals  0.5  (r*  =  CEP)  one 
must  integrate  equation  A. 10  with  respect  to  0  (0  to  2m)  and 

with  respect  to  r  (0  to  r*). 

This  leaves  : 

P(r  <  r*)  =  0.5  =  1  -  e xp [ -0 . 5 ( r */ s ) 2 1  (A. 11) 


and 


r  =  1.1774  S  =  CEP  ( A.  12  ) 

In  the  case  of  x  =  y  =  0  and  Sv  ^  Sw  the  formula  for 

x  y 


the  joint  distribution  is 


exp 


(r  cos9 \  2  /  r  sin9  \ 

— )  +  Hr) 


(A.  13) 


A- 5 


f ( r ,e) 


2  m  S  XS  y 
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If  0.25  <  Sx/Sy  <  1  this  expression  (after  integration) 

reduces  to: 

r*  =  0 . 6 1 4  S  y  +  0 . 5  6  3  S  y  =  CEP  (A. 14) 

(5:5).  For  other  values  of  Sx,  Sy  and  for  samples  not 
centered  on  the  target,  one  is  left  to  use  other 
approximations  to  determine  the  CEP. 

Popular  Approximation  Methods.  There  are  several 
approximation  methods  available  for  estimating  CEP  in 
elliptic,  biased  samples  (i.e.  Sx  ^  sy  and  x  4  y  4  0). 
Presented  here  are  four  of  the  approximation  methods 
available. 

G r ub bs -Pa tna i k/ Ch i -Sq ua re .  This  approximation  method 
is  based  on  the  work  of  Frank  E.  Grubbs  (2).  This  method 
uses  the  fact  that  the  bias  is  a  sum  of  noncentral  chi- 
squares  and  solves  for  the  radius  that  gives  a  probability 
of  C .  5  (i.e.  P[(x-x)2/Sx2  +  (y-y)2/Sy2  <=  R2 ]  =  0.5)  (2:55). 


This  is  the  approximation: 


CEP  = 


/ 


k v  /  2m 


(A.  15) 


where 


m  = 


k  = 


n  = 


( Sx2  +  Sy2  +  X2  +  y2 ) 

2(SX4  +  2p2Sx2Sy2  +  Sy4)  + 
4(x2Sx2  +  2xypSxSy  +  y2Sy2) 
F_1(0. 5) 

2  m  2  /  v 


A-  6 


(A.  16) 

(A.  17) 
(A.  18) 
(A. 19) 


and  F  is  the  chi-square  distribution  function  with  n  degrees 
of  freedom  (1:5).  In  this  formula  p  is  the  correlation 


coefficient.  If  the  correlation  coefficient  is  non-zero, 
one  can  perform  a  rotation  of  axes  into  a  coordinate  system 
where  the  two  coordinates  of  miss  distance  are  uncorrelated 
(3:20-21).  Attachment  1  to  this  tutorial  describes  how  to 
find  the  correlation  coefficient  and  how  to  perform  this 
rotation  to  eliminate  the  correlation  (this  makes  p  =  0). 

Grubbs-Patnaik/Wilson-Hilferty.  This  approximation 
technique  was  developed  as  a  modification  of  the  Grubbs- 
Pa tna ik/chi-sq uare  method.  The  W i 1 s on -H i 1 f e r t y  method 
transforms  the  chi-square  to  approximate  normal  variables. 
This  method  does  not  use  the  chi-square  function  described 
in  the  previous  section.  Here  is  the  G r u bb s -Pa tna i k / W i 1 s on - 
Hilferty  method: 

m{l-(v/(9m2)]3}  (A. 20) 

where  m  and  v  are  as  defined  in  the  G r u b b s -P a tna  ik / c h i - 
square  method  (equations  A. 16  and  A. 17)  (D:6). 

Modified  RAND-234.  This  method  is  a  fit  of  a  cubic 
polynomial  to  a  table  of  CEP  values.  Pesapane  and  Irvine 
used  regression  analysis  to  "...derive  a  mathematical 
expression  for  circular  error  probable  (CEP)  which 
approximates  probability  tables  contained  in  RAND  Report  R- 
234”  ( N : 2 ) .  This  is  the  Modified  RAND-234  method: 

CEP  =  CEPmpi(  1 .0039  -  0.0528v  +  0.4786v2  -  0.079  3v3)  (A. 21) 


where 


CEPflpi  =  0.614Ss  +  0.563Sij  (A. 22) 


The  Modified  RAND-234  was  developed  under  the  boundary 
condition  that  Sg/s^  is  greater  than  0.25  and  that  v  is  less 
than  or  equal  to  2.2  (1:2-3,  4:2-5).  If  one  has  a 
correlated  sample  as  described  in  the  section  on  the  Grubbs- 
Pat  na  i  k  /  c  h  i  -squa  re  method,  the  same  rotation  can  be  used  to 
remove  the  correlation. 

"Exact"  Method.  The  "exact”  method,  despite  its  name, 
is  also  an  approximation  technique.  It  is  recognized  as  the 
most  accurate  approximation  method  available.  This  method 
is  a  numerical  integration  technique  used  to  estimate  the 
integral  of  the  bivariate  normal  distribution.  This  method 
uses  Taylor  series  and  Infinite  series  expansion  to  estimate 
CEP.  The  "exact"  method  is  computationally  rigorous  and 
requires  a  computer  for  its  calculations.  Since  it  is  so 
computationally  rigorous  it  is  often  not  feasible  to  use  the 


exact"  method  and  that  Is  why  one  must  resort  to  u  s  i  n  ?  less 


accurate  approximations.  However,  the  "exact"  method  is 
useful  for  determining  which  approximation  technique  is  the 
most  accurate. 

The  "exact"  method  takes  the  bivariate  normal  function 
in  polar  form: 


P(R)  = 


R  2i  _  . 

)  /*“■>  -j( 

0  0  L  ' 


f(rcos@  -  x )  2  (rsin©  -  y)2 

—  + 


rdrd©  ( A. 2  7 ) 


where 


A  = 


2mSxSy 


(A. 28) 


After  series  expansion  this  is  the  equation  for  the  "exact’ 
approximation  technique: 


R2  a  oo  (  k  +  1 )  !  (  j  +1 )  ! 

P(R)  =  0.5  =  —  D  l  l  xkYj 


k=  0  j  =  0 


(k+j+1 ) ! 


(A. 29) 


where 


D  * 


SxSy 


exp 


x2 


y2 


2  S  x  2  .  2  S  y  2 


(2k) 


xk 


( k+1 ) ! ( k ! ) 2 


-R  2 


8SX2 


k 

V 

L 


k  ! 


1=0  ( k-1  )  !  ( 21 )  ! 


-  2  x  2 


(A  .  30) 


(A.  31  ) 


(  2  j  )  ! 


( j+l>! (j ! )2 


-R  2 


8Sy2 


j  ! 


I  - 

i  =  0  (j-i)! ( 2  i  )  ! 


-2y2 


(A. 32) 


(6:  1-5) 


( 


A  -  9 


To  compute  this  approximation,  one  must  first  make  a 


guess  at  the  CEP  (it  is  recommended  to  use  one  of  the  fast 
algebraic  approximation  methods  tc  arrive  at  your  guess), 
and  then  iterate  to  arrive  at  an  answer  that  is  accurate 
enough  for  one's  purpose.  This  method  sometimes  requires  a 
great  deal  of  computer  time  to  arrive  at  an  answer.  But  it 
is  the  most  accurate  estimator  of  CEP. 

Sample  Calculations 

Now  that  one  has  the  approximation  methods,  the  next 
step  is  to  actually  calculate  CEPs.  The  following  is  an 
example  of  how  one  would  do  this  for  each  of  the  presented 
approximation  techniques. 

Table  A-l  contains  sample  crossrange  and  downrange  miss 
distances  for  10  sample  impacts  of  a  theoretical  new  warhead 
(Note:  this  table  uses  feet  as  an  illustration,  the 

methods,  however,  are  insensitive  to  units  as  long  as  one  is 
consistent  with  the  units  used).  To  calculate  the  CEP  using 

Table  A-l 

Sample  Impacts  for  a  Theoretical  ICBM  Warhead 


Impact 

C  ros  s  range  Miss 
Dis  tance  ( x ) 

Down-ange  Miss 
Distance  (x) 

Radial  Miss 
Distance 

1 

400  f  t 

569  ft 

695.5  ft 

2 

324 

429 

545.6 

3 

-116 

125 

170.5 

4 

50 

-2  1  4 

219.8 

5 

-6  3 

-126 

143.2 

6 

257 

302 

396.6 

7 

76 

-  156 

173.5 

8 

96 

1  58 

184.9 

9 

-30 

53 

60 . 9 

10 

1  5  5 

204 

256  .2 

Since  the  correlation  coefficient  is  non-zero,  one  must 


perform  the  rotation  discussed  in  attachment  1.  After  this 
is  performed,  these  are  the  values  to  be  substituted  for  x, 
Sx,  y  ,  and  Sy: 


3T  =  37.1 

y'  =  281.8 

S'  =  91.1 

S'  =  290.7 

p  =  0 


Using  the  above  values  the  results  of  the  approximations 
are: 


Grubbs-Patnaik/chi-square: 

CEP  = 

324 . 34 

Grubbs-Patnaik/Wilson-Hilferty: 

CEP  = 

325.18 

Modified  RAND-234: 

CEP  = 

343.76 

"Exact" : 

CEP  = 

315.75 

Compare  these  values  to  the  non-para  metric  estimate,  202.35, 
and  one  can  see  how  inaccurate  the  non-parametric  estimate 
is,  given  a  small  sample  size.  The  calculations  required 
for  the  "exact"  method  can  often  require  a  great  deal  of 
computer  time  (anywhere  from  2  minutes  to  2  hours,  and 
sometimes  it  does  not  converge  to  an  answer  at  all),  so  it 
is  useful  to  know  which  other  approximation  methods  are  most 
accurate.  In  the  above  case  the  G r u bb s -Pa tna i k / c h i - s q u a r e 
approximation  technique  is  the  best  (closest  to  the  "exact") 


and,  or,  the  average,  the  Grutbs-Patnaik /  chi-square  is  the 
most  accurate  approximation  method. 


Ellipticity 
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Figure  A-4. 

Best"  Method  Across  the  Bias/Ellipticity  Parameter  Space 
/  =  RAND-234 

*  =  Gr u b bs -Pa tna i k/c h i -sq ua r e 

-  *  G r ubb s -Pa tna  i k /W i 1 s on-H i 1 f e r t y 

X  =  Areas  where  the  "exact"  does  not  converge 


However,  there  are  times  when  the  other  approximation 
methods  are  more  accurate  than  the  Grubbs-Pat n aik/chi- 
square.  Which  approximation  technique  is  best  varies  with 
the  bias  and  ellipticity  of  the  sample  impacts  (bias  is  the 
square  root  of  the  sum  of  the  squares  of  the  means  of  the 
crossrange  and  downrange  miss  distances,  ellipticity  is  the 
ratio  of  the  crossrange  standard  deviation  to  the  downrange 
standard  deviation).  Figure  A-4  shows  a  large  portion  of 
the  bias/ellipticity  parameter  space  and  which  approximation 
techniques  are  "best"  over  the  space.  The  technique  one 
will  ultimately  use  will  depend  on  how  accurately  one  needs 
to  estimate  the  CEP  of  a  weapon  system. 

It  should  be  noted  that  in  the  regimes  where  the 
"exact"  method  does  not  converge  there  are  other  techniques 
that  use  numerical  integration  techniques  to  estimate  CEP 
that  will  converge.  One  such  method  is  called  the 
correlated  bivariate  normal.  This  technique  is  not  as 
efficient  as  the  "exact"  but  it  will  converge  for  thos 
highly  elliptical  cases  where  the  "exact"  does  not. 

Conclusion 

Circular  error  probable  is  not  the  only  method  of 
representing  a  weapon  system's  accuracy,  however  it  is  a 
very  important  measure  that  is  widely  used.  Likewise  the 
approximation  methods  shown  here  are  not  the  only  methods 
available  for  estimating  CEP.  This  report  introduces  CEP 
and  describes  how  it  is  estimated. 


Attachment  1  to  CEP  Tutorial 

Correlation 

For  correlated  downrange  and  crossrange  miss  distances 
the  correlation  coefficient  must  be  taken  into 
consideration.  This  is  the  formula  used  to  calculate  the 
correlation  coefficient: 


n 

I  xiyi  ■  nJy 

i  =  1 


(n  -  l)(SxSy) 


(  A  .  1  .  1  ) 


If  the  correlation  coefficient  is  other  than  zero,  one  can 
perform  a  rotation  of  axes  into  a' coordinate  system  where 
the  two  coordinates  of  miss  distance  are  u n c o r r e 1 a te d .  This 
is  the  procedure: 


9  =  0.5{tan~1[2pSxSy/(Sx2  -  S  y  2 )  ] } 


(A  .  1  .  2  ) 


(S')2  =  0 . 5 [ ( S  2  +  S  2)  +  (S  2  -  Su2)/cos(2@)  ]  (A. 1.3) 


y  '  ’  '  “ x  "y 

(S')2  »  0 . 5 [ ( S  x  2  +  Sv2)  -  (sx2  -  Sy2)/ cos(29)  ]  (A. 1.4) 

(A. 1.5) 
(A.  1  .6) 


y/  V  •  -x  l  \  ^  ’  wy 

x '  =  xcos9  +  ysin9 


=  -xsin©  +  y  c  o  s  9 


Then  S',  S',  x' ,  y'  are  substituted  for  Sx,  Sy,  x,  y  in  the 
formulae  and  the  CEP  is  calculated  as  if  there  were  no 
correlation  (i.e.  p  =  0)  (3:20-21). 

Astbury  showed  that  the  formulae  involving  a 
correlation  coefficient  are  Invariant  to  this  rotation.  In 


A  .  I  -  1 


other  words,  they  give  the  sai?e  answer  whether  one  uses 


correlated  numbers,  or  the  uncorrelated  numbers  obtained 
from  the  rotation  (1:3-6). 


Attachment  2  to  CEP  Tutorial 
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Appendix  B 


Correlation 

For  correlated  downrange  and  crossrange  miss  distances 
the  correlation  coefficient  must  be  taken  into 
consideration.  This  is  the  formula  used  to  calculate  the 
correlation: 

n 

l  xiyi  -  nxy 
i  =  1 

p  = -  (B.l) 

(n  -  1)(S  Sv> 


If  the  correlation  coefficient  is  other  than  zero,  one  can 
perform  a  rotation  of  axes  into  a  coordinate  system  where 
the  two  coordinates  of  miss  distance  are  uncorrelated.  This 


is  the  procedure : 

G  =  0.5{tan_1[2pSxSy/(Sx2  -  Sy2)]}  (B.2) 
(S')2  =  0.5[(SX2  +  Sy2)  +  (Sx2  -  Sy2)/cos(20) ]  ( B . 3 ) 
(S')2  =  0 . 5 ( ( S  x  2  +  Sy2)  -  (Sx2  -  Sy2)/cos(26)  ]  (B.4) 
x '  =  xcos©  +  ysin©  (B.5) 
y '  =  -xsin©  +  ycos©  (B.6) 


Then  S't  sy,  x' ,  y'  are  substituted  for  Sx,  Sy,  x,  y  in  the 
formulae  and  the  CEP  Is  calculated  as  if  there  is  no 
correlation  (i.e.  p  =  0)  (3:20-21). 


B-l 


Astbury  showed  that  the  formulae  Involving  a 
correlation  coefficient  are  invariant  to  this  rotation.  In 
other  words,  they  give  the  same  answer  whether  one  uses  the 
correlated  numbers,  or  the  uncorrelated  numbers  obtained 
from  the  rotation.  (1:3-6). 


Appendix  C 


Secant.  Method 

The  Secant  method  is  used  to  approximate  the  root  of  a 
function.  In  this  case  the  root  is  the  CEP  (R),  and  the 
function  is  the  probability  ( P ( R ) ) .  The  secant  method  takes 
R^  and  Ri-i  (approximations  of  the  CEP)  and  takes  the  next 
approximation  (R^+j,  P(Ri+i))  as  the  intersection  of  the 
chord  joining  P(^i-l))  an<^  (&i»  F(Rj))  with  the  line 

y  =  0.5  (x  =  R).  This  is  the  formula  used: 

Ri  +  1  =  Rt  +  (Ri-i  -  R  i  )  (  P  i  + 1  -  Pi  >/ <Pi-l  "Pi)  (C-D 

where  P(R^  +  ^)  =  0.5  (1:9).  Figure  C-l  is  a  graphical 
representation  of  the  secant  method. 


Figure  C-l.  The  Secant  Method 

This  method  is  used  to  iterate  both  the  correlated 
bivariate  normal  and  the  "exact"  approximation  techniques  in 


this  study. 


Appendix  D 


Sample  Data 

This  appendix  contains  a  sample  of  the  data  collected 
for  this  study.  Table  D-l  contains  the  sample  cases  listing 
the  ellipticity,  bias,  sample  standard  deviations  and  sample 
means.  Table  D-2  is  a  listing  of  sample  CEPs  calculated  for 
each  of  the  approximation  techniques  for  each  case  listed  in 
table  D-l.  Finally  table  D-3  is  the  sample  relative  errors 
for  the  algebraic  approximations  and  the  CBN  relative  to  the 


exact"  method. 


CASE 

ELLIP 

S  IG  .  X 

S  IG  .  Y 

BIAS 

X  BAR 

Y  BAR 

1 

1  .00 

1000.00 

1000.00 

1555.63 

1100.00 

1100.00 

2 

1  .00 

1000 . 00 

1000 . 00 

14  14.21 

1000 .00 

1000 .00 

3 

1  .  00 

1000.00 

1000 . 00 

1  2  7  2.  7  9 

900 . 00 

900.00 

4 

1  .  00 

1000.00 

1000 . 00 

1131.37 

SCO  .  00 

800 .00 

5 

1.00 

1000 .00 

1000.00 

989  .  9  5 

700.00 

700.00 

6 

1  .  00 

1000.00 

1000 . 00 

848  .  5  3 

600 . 00 

600 .00 

7 

1  .00 

1000 .00 

1000 .00 

707  .  1  1 

500.00 

500 . 00 

8 

1 . 00 

1000 . 00 

1000.00 

565.69 

400.00 

400 .00 

9 

1  .00 

1000.00 

1000 .00 

424.26 

300.00 

300.00 

10 

1  .  00 

1000.00 

1000 . 00 

282.84 

200 .00 

200 .00 

1  1 

1 .00 

1000.00 

1000 . 00 

141.42 

ICO .00 

100 . 00 

12 

1 . 00 

1000 .00 

1000 .00 

1.41 

1  .  00 

1  .00 

1  3 

0 .90 

900.00 

1000.00 

1  5  5  5.6  3 

1100.00 

1100.00 

1  4 

0 .90 

900.00 

1000.00 

14  14.21 

1000 .00 

1000 .00 

15 

0 .90 

900.00 

1000 .00 

1272.79 

900 . 00 

900 . 00 

16 

0.90 

900 . 00 

1000 .00 

1131  .37 

800 . 00 

800 .00 

17 

0.90 

900.00 

1000 .00 

989  .  9  5 

700.00 

700.00 

13 

0.90 

900.00 

1000 . 00 

848  .  5  3 

600 . 00 

600.00 

19 

0.90 

900 .00 

1000.00 

707  .  1  1 

500.00 

500.00 

20 

0.90 

900 . 00 

1000 . 00 

565.69 

400.00 

400.00 

2  1 

0.90 

900 .00 

1000 . 00 

424.26 

300.00 

300.00 

22 

0 .90 

900.00 

1000 . 00 

282 .84 

200.00 

200.00 

23 

0.90 

900.00 

1000.00 

141.42 

100.00 

100 .00 

24 

0 .90 

900 . 00 

1000 . 00 

1.41 

1 . 00 

1  .00 

25 

0.80 

800.00 

1000 . 00 

155  5.6  3 

1100 .00 

1100 . 00 

26 

0.80 

800 . 00 

1000.00 

1414.2  1 

1000.00 

1000.00 

27 

0.80 

800.00 

1000.00 

1272.79 

900.00 

900.00 

28 

0.80 

800 . 00 

1000.00 

1131.37 

800 . 00 

800.00 

29 

0.80 

800.00 

1000.00 

989.95 

700.00 

700.00 

30 

0.80 

800.00 

1000 . 00 

848 .53 

600 . 00 

600 .00 

31 

0.80 

800 .00 

1000.00 

707  .  1  1 

500.00 

500.00 

32 

0.80 

800.00 

1000 . 00 

565.69 

400.00 

400.00 

33 

0 .80 

800.00 

1000 .00 

424.26 

300.00 

300 .00 

34 

0.80 

800.00 

1000 . 00 

282  .  84 

200.00 

200.00 

35 

0.80 

800.00 

1000.00 

141.42 

100.00 

100.00 

36 

0.80 

800 . 00 

1000 . 00 

1.41 

1  .  00 

1  .00 

37 

0.70 

700.00 

1000 .00 

1  5  5  5  .  6  3 

1100 . 00 

1100.00 

38 

0.70 

700.00 

1000.00 

14  14.21 

1000 . 00 

1000 .00 

39 

0.70 

700.00 

1000 .00 

1272.79 

900.00 

900.00 

40 

0 .70 

700 . 00 

1000.00 

1131.37 

800.00 

800 .00 

41 

0.70 

700 . 00 

1000.00 

989.95 

700 . 00 

700 .00 

42 

0 .70 

700 .00 

1000 . 00 

848  .  53 

600 . 00 

600 .00 

43 

0 .70 

700.00 

1000 .00 

707  .  1  1 

500 . 00 

500.00 

44 

0.70 

700 . 00 

1000 . 00 

565.69 

400.00 

400.00 

45 

0 .70 

700.00 

1000.00 

424  .  26 

300.00 

300 .00 

46 

0.70 

700.00 

1000 .00 

282.84 

200 . 00 

200.00 

47 

0.70 

700.00 

1000.00 

141.42 

100 .00 

100.00 

48 

0. 70 

700.00 

1000.00 

1.41 

1  .  00 

1  .00 

Table  D-l 

Sample  Cases  (part  1) 


D-2 


CASE 

ELLIP 

SIC.  X 

S  IG  .  Y 

BIAS 

X  PAR 

Y  BAR 

49 

0.60 

600.00 

1000 .00 

1555.63 

1100.00 

1100.00 

50 

0 .60 

600 . 00 

1000.00 

1414.21 

1000 . 00 

1000 .00 

51 

0 .60 

600.00 

1000 . 00 

1272.79 

900 .00 

900 . 00 

52 

0 .60 

600 . 00 

1000 . 00 

1131  .37 

800 .00 

800.00 

53 

0.60 

600 .00 

1000.00 

989.95 

700 .00 

700.00 

54 

0.60 

600.00 

1000 . 00 

848.53 

6C0 . 00 

600 . 00 

55 

0.60 

600.00 

1000.00 

707.11 

500.00 

500.00 

56 

0.60 

600.00 

1000.00 

56  5.69 

400.00 

400.00 

57 

0.60 

600.00 

1000 . 00 

424.26 

300 .00 

300 . 00 

58 

0 .60 

600 . 00 

1000 .00 

282  .  84 

200.00 

200.00 

59 

0.60 

600.00 

1000.00 

141.42 

100.00 

100 . 00 

60 

0 .60 

600 . 00 

1000 . 00 

1.41 

1  .  00 

1  .00 

6  1 

0 . 50 

500 .00 

1000 .00 

1  5  5  5  .  6  3 

1100.00 

1100.00 

62 

0 . 50 

500.00 

1000 . 00 

14  14.21 

1000 . 00 

1000 .00 

63 

0 .50 

500.00 

1000 . 00 

1272.79 

900.00 

900.00 

64 

0 . 50 

500.00 

1000 . 00 

1131.37 

800 . 00 

800.00 

65 

0 .50 

500.00 

1000 . 00 

989.95 

700.00 

700.00 

66 

0 . 50 

500.00 

1000 . 00 

848  .  5  3 

600 . 00 

600 . 00 

67 

0.50 

500 .00 

1000 . 00 

707  .  1  1 

500 . 00 

500 . 00 

68 

0 . 50 

500.00 

1000 . 00 

565.69 

400 .00 

400.00 

69 

0.50 

500.00 

1000.00 

424.26 

300.00 

300 . 00 

70 

0.50 

500.00 

1000 . 00 

282 .84 

200.00 

200 .00 

71 

0 .50 

500.00 

1000.00 

141.42 

100.00 

100.00 

72 

0 . 50 

500.00 

1000.00 

1.41 

1.00 

1 .00 

73 

0.40 

400 . 00 

1000.00 

1  555.6  3 

1100 . 00 

1100 . 00 

74 

0.40 

400.00 

1000 . 00 

1414.21 

1000 .00 

1000.00 

75 

0.40 

400.00 

1000.00 

1272.79 

900.00 

900.00 

76 

0.40 

400 . 00 

1000 . 00 

1131.37 

800 . 00 

800.00 

77 

0.40 

400 .00 

1000.00 

989.95 

700.00 

700.00 

78 

0 . 40 

400 . 00 

1000 . 00 

848.53 

600.00 

600.00 

79 

0.40 

400 .00 

1000.00 

707  .  1  1 

500.00 

500.00 

80 

0.40 

400.00 

1000 .00 

565.69 

400 .00 

400.00 

8  1 

0.40 

400 .00 

1000 . 00 

424  .  2  5 

300 . 00 

300 . 00 

82 

0 . 40 

400.00 

1000 . 00 

282  .  84 

200.00 

200.00 

83 

0.40 

400.00 

1000.00 

141.42 

100.00 

100 . 00 

84 

0.40 

400 . 00 

1000 . 00 

1.41 

1  .00 

1 .00 

85 

0 . 30 

300 . 00 

1000.00 

1555.63 

1100.00 

1100.00 

86 

0.30 

300.00 

1000 . 00 

1414.21 

1000.00 

1000 . 00 

87 

0 .30 

300 .00 

1000 . 00 

1272.79 

900 .00 

900.00 

88 

0 . 30 

300 . 00 

1000 . 00 

1131.37 

800.00 

800 .00 

89 

0.30 

300.00 

1000 . 00 

989  .  9  5 

700.00 

700 . 00 

90 

0.30 

300 . 00 

1000 . 00 

848.53 

600.00 

600 . 00 

9  1 

0 .30 

300.00 

1000.00 

707  .  1  1 

500 .00 

500 . 00 

92 

0 . 30 

300 . 00 

1000 . 00 

565.69 

400 . 00 

400 . 00 

93 

0.30 

300.00 

1000.00 

424.26 

300.00 

300.00 

94 

0 . 30 

300 . 00 

1000 . 00 

282  .  84 

200.00 

200.00 

95 

0.30 

300.00 

1000.00 

141.42 

100 .00 

100.00 

96 

0 . 30 

300 . 00 

1000.00 

1.41 

1  .00 

1  .00 

Table  D-I 

Sample  Cases  (pari  2) 


I? 


20 

20 

20 

20 

20 

20 

20 


20 

20 

20 


10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

10 

01 

01 

01 

01 

01 

01 

01 

01 

01 

01 

01 

01 


S  IG  .  X 

200 .00 
200.00 
200.00 
200.00 
200 .00 
200 . 00 
200.00 
200 .00 
200.00 
200 . 00 
200 .00 
200.00 
100 . 00 
100.00 
100 .00 
100.00 
100.00 
100.00 
100.00 
100.00 
100 .00 
100 . 00 
100 .00 
100 . 00 
10 . 00 
10.00 
10.00 
10.00 
10.00 
10.00 
10 .00 
10.00 
10.00 
10.00 
10.00 
10 .00 


SIG  .  Y 

1000 .00 
1000 . 00 
1000 .00 
1000 . 00 
1000 .00 
1000 . 00 
1000.00 
1000 .00 
1000 .00 
1000 . 00 
1000 .00 
1000 . 00 
1000.00 
1000 . 00 
1000 .00 
1000 . 00 
1000.00 
1000 . 00 
1000.00 
1000.00 
1000 .00 
1000 . 00 
1000.00 
1000 . 00 
1000.00 
1000.00 
1000 . 00 
1000 . 00 
1000.00 
1000 . 00 
1000.00 
1000 . 00 
1000 . 00 
1000 . 00 
1000 . 00 
1000.00 


BIAS 

1555.63 

1414.21 

1272.79 

1131.37 
989.9  5 

848.53 
707  .  1  1 

565.69 
424.26 
282.84 

141.42 

1.41 
1555.63 

1414.21 

1272.79 

1131.37 
989  .  9  5 

848.53 
707  .  1  1 

565.69 
424.2  6 
282  .84 

141.42 

1.41 
1555.6  3 

1414.21 

1272.79 

1131.37 
989.9  5 

848.53 
707  .  1  1 

565.69 
424.26 
282.84 

141.42 

1.41 


X  BAR 

1130.00 
1000 .00 
900 .00 
800.00 
700.00 
600 . 00 
500.00 
400 .00 
300 . 00 
200 . 00 
100.00 
1 . 00 
1100.00 
1000 . 00 
900.00 
800.00 
700 . 00 
600 . 00 
500.00 
400 . 00 
300.00 
200 . 00 
100 . 00 
1.00 
1100.00 
1000.00 
900.00 
800 . 00 
700 . 00 
600 . 00 
500 .00 
400 . 00 
300.00 
200 .00 
100 .00 
1 . 00 


Y  EAR 

1100.00 
1000 .00 
900 .00 
800.00 
700 . 00 
600 .00 
500.00 
400.00 
300 . 00 
200 .00 
100 . 00 
1  .  00 
1100.00 
1300.00 
900 . 00 
800 .00 
700.00 
600.00 
500.00 
400 .00 
300.00 
200 .00 
100.00 
1.00 
1100.00 
1000 . 00 
900.00 
800 . 00 
700.00 
600 . 00 
500.00 
400.00 
300.00 
200.00 
100 . 00 
1  .  00 


■  Table  D-l 
Sample  Cases  (part  3) 


D-4 


CASE 

Grubbs/ CS 

Grubbs/ WH 

1 

1857.29 

1861  .86 

2 

1750.71 

1755.28 

3 

1650.01 

1654.75 

4 

1556.27 

1561.32 

5 

1470.70 

1476.14 

6 

1  394 .53 

1400.35 

7 

1328.91 

1335.03 

8 

1274.75 

1281 .04 

9 

1232.61 

1238.95 

10 

1202 .70 

1208.99 

1  1 

1184.89 

1191.12 

12 

1178.98 

1185.19 

13 

1825.38 

1829.81 

14 

1715.94 

1720.28 

15 

1612.11 

1616.56 

16 

1515.04 

1519.76 

17 

1426.00 

1431.10 

18 

1  346  .  36 

1351.88 

19 

1277.44 

1283.31 

20 

1220 . 34 

1226.45 

2  1 

1175.84 

1182.01 

22 

1144.46 

1150.34 

23 

1125.69 

1131.43 

24 

1119.45 

1125.16 

25 

1795.71 

1799.82 

26 

1683  .  10 

1687.27 

27 

1576.01 

1580.23 

28 

1475.35 

1479.84 

29 

1382.49 

1387.41 

30 

1298.90 

1304.33 

31 

1226  .  1  1 

1232.01 

32 

1165.71 

1171  .72 

33 

1118.96 

1124.38 

34 

1085.45 

1090.52 

35 

10 65.37 

1070.28 

36 

1058.70 

1063.55 

37 

1768  .4  1 

1772.17 

38 

1652.56 

1656.61 

39 

1  542  .  1  1 

1546.17 

40 

1437.69 

1442.04 

41 

1340.69 

1345.55 

42 

1252.70 

1258.23 

43 

1175.51 

1181.65 

44 

1112.22 

1117.32 

45 

1062.04 

1066.47 

46 

1025.85 

1029.88 

47 

1004 . 08 

1007.90 

48 

996.82 

1000 . 59 

Table 
Sample  CEP 


RAND  CBN  Exact 

1867.99  1850.30  1869.75 

1758.27  1740.69  1 652.25 

1655.09  1637.23  1654.24 

1559 .44  1541.35  1557.31 

1472.28  1454.52  1469.56 

1394.59  1378.04  1392.28 

1327.33  1312.84  1326.38 

1271.48  1259.46  1272.47 

1228.01  1218.10  1230.69 

1197.89  1188.72  1201.00 

1182.09  1171.18  1183.28 

1x81.52  1165.35  1177.38 

1836.14  1820.39  1839.29 

1723.07  1707.95  1669.30 

1616.36  1601.20  1555.53 

1517.07  1501.66  1517.38 

1426.29  1410.96  1425.86 

1345.11  1330.64  1344.82 

1274.59  1261.90  1275.46 

1215.83  1205.50  1218.55 

1169.91  1161.75  1174.39 

1  1  37.90  1  130.67  1  143. CO 

1120.88  1112.11  1124.27 

1119.88  1105.94  1118.04 

1806.21  1791.28  1809.55 

1689.80  1676.48  1643.36 

1579.45  1566.84  1583.02 

1476.35  1463.87  1416.66 

1381.73  '  1369.27  1383.93 

1296.79  1234.81  1298.86 

1222.74  1212.02  1225.53 

1160.80  1151.97  1165.03 

1112.18  1105.25  1117.93 

1078.08  1071.99  1084.38 

1059.72  1052.12  1064.34 

1058.24  1045.51  1057.67 

1778.26  1762.20  1779.77 

1658.58  1645.52  1618.01 

1544.53  1533.56  1549.18 

1437.48  1427.66  1442.56 

1338.78  1329.46  1280.21 

1249.81  1240.81  1254.65 

1  1  7  1  .  9  3  1  16  3.56  1  1  7  6.99 

1106.50  1099.23  1112.28 

1054.90  1048.80  1061.52 

1018.48  1012.74  1025.20 

998.62  991.14  1003.44 

996.60  983.95  996.19 

D-2 

s  (part  1  ) 


E Si 


D-5 


rubbs/ CS 

Grubbs/WH 

RAND 

CBN 

Exact 

1743.70 

1747.20 

1752.28 

1732.18 

1749.07 

1624  .  7  3 

1628.68 

1629.51 

1614.08 

1592.31 

1510.93 

1514.88 

1511.78 

1500.40 

1515.  39 

1402.66 

1406.94 

1400 .65 

1392.25 

1  406.60 

1301.30 

1306.23 

1297.65 

1291  .09 

1  304.9  3 

1208.51 

1214.36 

1204.36 

1198.64 

1146.17 

1127.88 

1133.08 

1122.30 

1116.87 

1130.13 

1060 .05 

1064.15 

1053.05 

1047.72 

1060.74 

1005.74 

1009 . 14 

998.14 

99  2.7  9 

1005 . 58 

966.23 

969.22 

959.14 

9  53  .  1  3 

96  5.7  1 

942.30 

945.09 

937.58 

929  .  2  2 

941.67 

9  34.30 

937.03 

934.96 

921.25 

9  33  .  6  4 

1721.93 

1725.24 

1728.11 

1700.22 

1716.59 

1600 . 06 

1603.94 

1602.62 

1531.01 

1  56  5  .  1  8 

1483.04 

1486 . 9  1 

1431.35 

1466.07 

1480.50 

1370.97 

1375.24 

1366.03 

1356.35 

1  370  .  16 

126  5  .  1  9 

1270.30 

1258.59 

1253.05 

1  2  66.4  5 

1167.8  1 

1173.76 

1160.67 

1157  .63 

1170.78 

1083.23 

1087.50 

1074.07 

1071.87 

1013.89 

1010.35 

1013.54 

1000 .60 

997.88 

1010.84 

951.25 

9  5  3.8  5 

942.01 

9  37  .  85 

950.74 

907.7  5 

910.08 

900  .  10 

89  3.70 

906 .50 

881.15 

88  3.40 

876.63 

366.7  5 

879.48 

872.2  1 

874.44 

873.32 

857 .70 

870.40 

1703.48 

1706.65 

1705 . 39 

1665.27 

17  0  6.6  5 

1579.03 

1582.85 

1577.81 

1545.09 

1  5  3  5  .  1  4 

1459 .04 

1462.86 

1453.33 

1429  .  19 

1443.34 

1343.4  1 

1347.72 

1333.99 

1318.40 

1331.93 

1233 . 38 

1238.78 

1221.87 

1213.71 

1226.77 

1132.56 

1137.69 

1119.01 

1116.30 

1129.14 

1042.99 

1046.46 

102  7  .  48 

1027.65 

1040.45 

964.77 

967.35 

949 .33 

949 .67 

881.86 

900.40 

902 .7  3 

886.63 

884.7  4 

397.83 

8  52  .  33 

854.79 

841.43 

835.62 

348.32 

822.60 

825.29 

815.78 

804 .92 

8  18.21 

812.54 

815.34 

811.68 

794.47 

807.79 

1688.69 

1691.76 

1683.33 

1626.27 

1665.70 

1562 . 08 

1565.87 

1554.75 

1505.00 

1544.35 

1439 .57 

1443.36 

1427.69 

1388.25 

1437.18 

1  320.83 

1325.22 

1304 . 56 

1388.25 

1307.55 

1207.01 

1212.01 

1187.76 

1172.30 

1205.28 

1103.13 

1107.56 

1079 . 7  1 

1075.52 

1085.98 

1008 .86 

1011.76 

982.82 

973.40 

994  .  33 

925.43 

9  27.8  1 

899.50 

888.73 

9  15.50 

855.7  3 

858.44 

832.16 

821.45 

74  5.28 

802.87 

806.41 

783.21 

774.59 

794.54 

769.7  3 

774.10 

755.06 

740.1  3 

760.92 

758.43 

763.14 

750.04 

738.92 

749.96 

Table  D-2 

Sample  CEPs  (part  2) 


rubbs / CS 

Grubbs/ WH 

RAND 

CBN 

Exact 

1677.88 

1630.89 

1660.44 

1585.93 

* 

1549.65 

1553.40 

1532.65 

1462.39 

* 

1425.18 

1428.96 

1404. 17 

1344 . 02 

* 

1304.02 

1308.49 

1277.84 

1231  .65 

* 

1187.17 

1193.22 

1156.54 

1125.89 

* 

1080.84 

1084 .78 

1043. 13 

1027.02 

* 

982.67 

985.27 

940.46 

936.56 

* 

894 . 78 

897.30 

851.40 

856  .  12 

870.61 

820.35 

823.93 

778.81 

788.00 

802 . 8  3 

763.08 

768.39 

725.56 

7  3  5.30 

607 . 9  5 

726.72 

733.65 

694.49 

701.57 

7  17.34 

714.23 

721.81 

688.40 

689.91 

705.77 

1671.29 

1674.26 

1633.90 

1  560.2  2 

1653.24 

1542.04 

1545.78 

1509.86 

1433.68 

it 

141o.35 

1420.12 

1381.96 

1312.60 

* 

1293.63 

1298.17 

1253.66 

1198.19 

* 

1175.02 

1181.12 

1128.40 

1091.35 

* 

1066.91 

1070.56 

1009.65 

992.36 

it 

966  .  1  3 

968 . 6  1 

900.84 

902 .83 

* 

875.20 

877.98 

805.43 

824  .  15 

* 

797.44 

801.93 

726.87 

7  5  8.36 

* 

736.99 

744.01 

668.62 

708  .  1  3 

it 

698.26 

707.60 

634.12 

676.38 

411.30 

684.88 

695.16 

626.76 

665.49 

615.00 

1669.09 

1672.06 

1602.61 

1667.58 

* 

1539.51 

1543.24 

1486 . 09 

1572 . 47 

* 

1413.39 

1417.17 

1361.66 

1402 .40 

* 

1290.15 

1294.72 

1233.47 

1319.39 

* 

1171.07 

1177.05 

1105.68 

1192.82 

* 

1062.22 

1065.77 

982.45 

1066 . 58 

* 

960.53 

962.98 

367.93 

977.34 

it 

868  .5  1 

871.43 

766.27 

863.73 

* 

789 . 57 

794.43 

681.63 

751.08 

* 

727.97 

735.68 

618.16 

817.96 

* 

688.36 

698.68 

580.02 

774.81 

it 

674.64 

686 . 03 

571.29 

570.98 

it 

Table  D-2 

Sample  CEPs  (part  3) 

Note:  *  =  "exact"  does  not  converge 


L-  W  V"  „T» 

r.'  v.  w  -  .•  v 

CASE 

Grubhs/CS  G 

rubbs /WH 

RAND 

cp:; 

i 

"PEST"  £■ 

1 

-0 . 0067 

-0 .0042 

-0 . 0009 

-0.0104 

2 

0.0596 

0.0624 

0.0642 

0.05  3  5 

1 

3 

-0 . 0026 

0 .0003 

0 . 0005 

-0.0103 

2  fi 

4 

-0.0007 

0 . 0026 

0.0014 

-0.0102 

1  P 

5 

0 . 0008 

0 . 0045 

0.0019 

-0.0102 

i  •: 

6 

0.0016 

0 .0058 

0.0017 

-0.0102 

i  :■ 

7 

0.0019 

0 . 0065 

0 .0007 

-0.0102 

3  ;; 

8 

0.0018 

0.0067 

-0 . 0008 

-0.0102 

3  ►; 

9 

0.0016 

0 . 0067 

-0.0022 

-0.0102 

1  i 

10 

0.0014 

0.0067 

-0 . 0026 

-0.0102 

i  p 

1  1 

0.0014 

0 .0066 

-0.0010 

-0.0102 

3  *: 

12 

0.0014 

0 . 0066 

0.0035 

-0.0102 

1 

1  13 

-0 . 0076 

-0 . 0052 

-0.0017 

-0 .0103 

3 

,  14 

0.0279 

0.0305 

0.0322 

0.0232 

i  •; 

15 

0.0364 

0.0392 

0.039  1 

0 .0294 

16 

-0.0015 

0.0016 

-0.0002 

-0.0104 

3  P 

1  7 

0.0001 

0 . 0037 

0 . 0003 

-0.0104 

l  •! 

:  18 

0.0011 

0.0053 

0 . 0002 

-0.0105 

3 

!  19 

0.0015 

0 .0062 

-0 . 0007 

-0 .0106 

3 

1  20 

0.0015 

0.006  5 

-0.0022 

-0.0107 

i  •; 

i  21 

0.0012 

0.0065 

-0 . 0038 

-0.0108 

1  i 

1  22 

0.0013 

0.0064 

-0.0045 

-0.0108 

i  f 

!  23 

0.0013 

0 .0064 

-0 .0030 

-0.0108 

i 

'  24 

0.0013 

0.0064 

0.0016 

-0.010 8 

l  V 

!  25 

-0 . 0076 

-0 . 0054 

-0.0018 

-0.0101 

3 

26 

0 .0242 

0.0267 

0.0283 

0 .0202 

1 

27 

-0 . 0044 

-0 . 0018 

-0.0023 

-0 .0102 

2  1 

28 

0.0414 

0.0446 

0.042  1 

0.0333 

1  P 

'!  2  9 

-0 . 0010 

0 . 0025 

-0.0016 

-0.0106 

1 

!  30 

0.0000 

0.0042 

-0.0016 

-0.0108 

1 

;  3 1 

0 .0005 

0 . 0053 

-0 . 0023 

-0.0110 

1  s’ 

i  32 

0.0006 

0.0057 

-0.0036 

-0.0112 

1  £ 

33 

0 . 0009 

0 . 0058 

-0.0051 

-0.0113 

1  1 

34 

0.0010 

0.0057 

-0 . 0058 

-0.0114 

1  " 

35 

0.0010 

0.0056 

-0.0043 

-0.0115 

1  ^ 

;  36 

0.0010 

0.0056 

0 . 0005 

-0.0115 

3 

1  37 

-0.0064 

-0 .0043 

-0 . 0008 

-0 . 0099 

3  N 

38 

0.0214 

0.0239 

0.025  1 

0.0170 

1  N 

39 

-0 .0046 

-0.0019 

-0 . 0030 

-0.0101 

2  1 

40 

-0.0034 

-0.0004 

-0.0035 

-0.C103 

2 

41 

0.0472 

0.0510 

0.0457 

0.0385 

3 

42 

-0.0016 

0.0028 

-0 . 0039 

-0.0110 

i 

43 

-0.0013 

0.0040 

-0 . 0043 

-0.0114 

i 

44 

-0.0000 

0.0045 

-0 . 0052 

-0.0117 

1  N 

4  5 

0.0005 

0 . 0047 

-0 . 0062 

-0.0120 

1  I 

46 

0.0006 

0 .0046 

-0.0066 

-0.0122 

1  -3 

47 

0 . 0006 

0 . 0044 

-0 . 0048 

-0.0123 

i 

48 

0.0006 

0.0044 

0.0004 

-0.0123 

3 
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CASE 

Grubbs / CS 

Grubbs /WH 

RAND 

CBN 

"BEST 

49 

-0.0031 

-0.0C11 

0.0018 

-0 .0097 

2 

50 

0 .0204 

0.0229 

0.0234 

0.0137 

1 

51 

-0 . 0029 

-0 . 0003 

-0 . 0024 

-0 . 0099 

2 

52 

-0 . 0028 

0.0002 

-0.0043 

-0.0102 

2 

53 

-0 .0028 

0.0010 

-0.0056 

-0.0106 

2 

54 

0 .0544 

0.0595 

0.0508 

0 . 0458 

3 

55 

-0 . 0020 

0 . 0026 

-0 . 0069 

-0.0117 

1 

56 

-0 .0006 

0.0032 

-0.007  2 

-0.0123 

1 

57 

0 .0002 

0.0036 

-0 . 007  3 

-0.0126 

1 

58 

0 . 0006 

0.0037 

-0.0068 

-0.0130 

1 

59 

0 . 0007 

0.0036 

-0.0043 

-0.0132 

1 

60 

0 . 0007 

0 . 0036 

0.0014 

-0.0133 

1 

61 

0.003  1 

0.0050 

0 . 0067 

-0.0095 

1 

62 

0.0223 

0 .0248 

0.0239 

0.0101 

1 

63 

0.0017 

0 . 0043 

0 .0006 

-0 .0097 

3 

64 

0 . 0006 

0.0037 

-0 . 0030 

-0.0101 

1 

65 

-0.0010 

0.0030 

-0 . 0062 

-0.0106 

1 

6  6 

-0 . 0025 

0.0025 

-0.0086 

-0.0112 

1 

67 

0 .0684 

0.0726 

0.0594 

0.0572 

3 

68 

-0.0005 

0.0027 

-0.0101 

-0 .0128 

1 

69 

0 . 0005 

0.0033 

-0 . 0092 

-0.0136 

1 

70 

0.0014 

0 . 0040 

-0.0071 

-0.0141 

1 

7  1 

0.0019 

0.0045 

-0 . 0032 

-0.0145 

1 

72 

0.0021 

0.0046 

0.0034 

-0.0146 

1 

73 

-0 . 0019 

0.0000 

-0 . 0007 

-0.02^2 

2 

74 

0.0286 

0.0311 

0.0278 

0.0065 

3 

75 

0.0109 

0.0135 

0 . 0069 

-0.0098 

3 

76 

0.0086 

0.0119 

0.0016 

-0.0102 

3 

77 

0.0054 

0 . 0098 

-0 . 0040 

-0 .0106 

3 

78 

0 . 0030 

0.0076 

-0.0090 

-0.0114 

1 

79 

0 .0024 

0 . 0058 

-0.0125 

-0.0123 

1 

80 

0.0940 

0 . 0969 

0.0765 

0.0769 

3 

81 

0.0029 

0.0055 

-0.0125 

-0 .0146 

1 

8  - 

0.0041 

0.0070 

-0.0087 

-0.0156 

1 

83 

0 . 0054 

0 . 0087 

-0 .0030 

-0 .0162 

3 

84 

0.0059 

0.0093 

0.0048 

-0.0165 

3 

85 

0.0138 

0.0156 

0.0106 

-0.0237 

3 

86 

0.0115 

0.0139 

0 . 0067 

-0.0255 

3 

87 

0.0017 

0.0043 

-0.0066 

-0 .0340 

1 

88 

0.0102 

0.0135 

-0.0023 

-1 .0000 

3 

89 

0.0014 

0.0056 

-0.0145 

-0 . 0027 

1 

90 

0.0158 

0.0199 

-0.0058 

-0 . 0096 

3 

9  1 

0.0146 

0.0175 

-0.0116 

-0 .02  10 

3 

92 

0.0108 

0.0134 

-0.0175 

-0.0292 

1 

93 

0.1482 

0.1518 

0.1166 

0 . 1022 

3 

94 

0.0105 

0.0149 

-0.0143 

-0.0251 

1 

95 

0.0116 

0.0173 

-0 .0077 

-0.0273 

3 

96 

0.0107 

0.0176 

0.0001 

-0.0147 

3 

Table  D-3 

Sample  Relative  Errors  (part  2) 


D-9 


Grubbs/CS 


Grubbs/WK 
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0.02  78 
0.0218 
0.2552 
0.0131 
0.0120 
0 . 0079 


0.0307 
0 .0263 
0.2639 
0 .0227 
0.0227 
0 . 0097 


-0.0221  -0.0166  3 

-0.0058  -0.0096  3 

0.1934  0.2095  3 

-0.0319  -0.0220  1 

-0.0246  -0.0225  1 

-0.0147  -0.0591  1 


0.6983  0.7202  0.5426 

0. 1  1  36  0.  1  303  0.019  1 


0.6448  3 

0.0821  3 
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Appendix  E 


This  appendix  contains  the  source  code  for  the  computer  program  developed 
by  this  study.  The  program  is  written  in  Turbo  Pascal  (Turbo  Pascal  is  a 
product  of  Eorland  International),  a  version  of  the  Pascal  computer  programming 
language  for  personal  computers.  The  program  was  written  and  run  on  a  SANYO 
MBC-550  personal  computer  (operating  system:  MS-DOS)  owned  by  the  author.  This 
program  was  also  run  on  a  Zenith  Z-150  personal  computer  at  the  Air  Force 
Institute  of  Technology.  With  some  minor  alterations  this  program  can  be  run  on 
any  computer  with  any  version  of  the  Pascal  programming  language. 

This  program  runs  interactively  with  the  user  inputs  to  calculate  circular 
error  probable  (CEP)  using  five  different  approximation  techniques.  The  program 
will  calculate  CEP  for  all  five  methods,  only  the  methods  the  user  selects,  or 
it  will  select  a  "best"  method  from  among  the  three  algebraic  methods  which  give 
results  in  seconds.  The  user  decides  which  method(s)  the  program  will  compute. 


E-l 


{THIS  PROGRAM  WAS  DEVELOPED  BY  CAPT  RICHARD  L.  ELDER  AS  A  PART  OF  A 


MASTER'S  THESIS  AT  THE  AIR  FORCE  INSTITUTE  OF  TECHNOLOGY. 

TITLE  OF  THE  THESIS  IS  "AN  EXAMINATION  OF  CIRCULAR  ERROR  PROBABLE 
APPROXIMATION  TECHNIQUES."} 


program  CEPAPPROVIMATIONANALYSIS  ; 

var  sigx,  sigy,  mux,  muy,  rho,  ellip,  bias:  real; 

cbnprob,  randc,  gpcs,  gpwh,  prob,  cep,  cbn:  real: 
pick,  ints:  integer; 


********** 


function  atotheb(a,b:  real):  real; 

{This  function  takes  an  input  variable  a  and  raises  it 
to  the  b  power} 


begin 

if  a  =  0  then  atotheb:=  0 

else  if  (b*ln(a))  <  -70  then  atotheb:3  0 
else  atotheb:=  exp(b*ln(a)) 
end  {atotheb}  ; 


********** 


function  interpolated ,  aup,  alow,  kup ,  klow:  real):  real; 

{This  function  is  a  function  to  perform  linear 

interpolation.  It  is  used  in  procedure  GPCSCEP 
to  interpolate  between  values  found  in  a  table 
of  the  chi-square  function.} 

var  percent:  real; 

begin 

percent:=  (a  -  alow)/(aup  -  alow); 
interpolate:3  klow  +  (percent*(kup  -  klow)) 
end  {interpolate}; 


********** 


procedure  GPCSCEP(var  sigx,  sipy,  mux,  muv,  rho:  real); 

{This  procedure  calculates  CEPs  according  to  the 

Grubbs-Patnaik/chi-square  approximation  technique) 

var  m ,  x,  y,  v,  n,  k:  real; 

begin 

m:=  sqr(sigx)+sqr(sigy)+sqr(mux)+sqr(muy) ; 
x:=  atotheb(sigx,  4); 
y:=  atotheb(sigy,  4); 
v:=  2*(x+(2*sqr(rho*sigx*sigy) )+y) ; 

v:=  v+(4*(sqr(mux*sigx)+(2*mux*muy*rho*sigx*sigy)+sqr(muy*sigy) ) ) ; 

n:=  2*sqr(m)/v; 

if  n  >=  11  then  k:=  (n-0.7); 

case  trunc(n)  of 

0:  k:=  interpolate(n, 1 ,0,0.455,0) ; 

1:  k:=  interpolate(n , 2 , 1  ,  1 . 39  ,0. 455) ; 

2:  k:=  interpolate (n , 3 , 2 , 2 . 37 , 1 . 39) ; 

3:  k:=  interpolate(n ,4 , 3 , 3 . 36 ,2 . 37 ) ; 

4:  k:=  interpolated , 5 , 4 ,4 . 35 , 3 . 36) ; 

5:  k:=  interpolate (n , 6 , 5 , 5 . 35 ,4 . 35 ) ; 

6:  k:=  interpolate(n , 7 , 6 , 6 . 35 , 5 . 35) ; 

7:  k:=  interpolated ,  8 , 7 , 7 . 34 , 6 . 35) ; 

8:  k:=  interpolated  ,  9 , 8 , 8 . 34 ,7 . 34) ; 

9;  k:=  interpolated ,  10, 9 ,9 . 34 ,8 . 34) ; 

10:  k:=  interpolated,  11 , 10, 10. 3,9. 34) 
end ; 

gocs:  =  sqrt( (k*v)/(2*m) ) 
end  {GPCSCEP}; 


********** 


procedure  GPWHCEP(var  sigx,  sigy,  mux,  muy:  real); 

(This  procedure  calculates  CEPs  according  to  the 
Grubbs-Patnaik/Wilson-Hilf ertv  approximation 
technique. ) 

var  aa,  m,  v,  vv,  wv:  real; 
begin 

m:=  sqr(sigx)  +  sqr(sigy)  +  sqr(mux)  +  sqr(muy); 

vv:=  sqr(sqr(sigx) )  +  2*  rho  *  sqr(sigx)  *  sqr(sigy)  +  sqr ( sq r ( s igy ) ) ; 
vvv:=  sqr(mux)*sqr(sigx)  +  2*mux*muy*rho*sigx*sigv  +  sq r(muy )*sq r ( s 1 gy ) ; 
v:=  2*vv  +  4*vvv; 
aa:=  1  -  (v/(9*sqr(m) ) ) ; 
gpwh:=  sqrt(m*(atotheh(aa,  3))); 
end  {GPWHCEP} ; 


•kirkic’k’k'k'k'k'k 


procedure  RANDCEP  (var  sigx,  sigy,  mux,  muy:  real); 


{This  procedure  calculates  CEPs  according  to  the 
modified  RAND-234  approximation  technique} 


var  radical,  sigxy,  sigS,  sigL,  b,  v,  cepmu:  real; 
begin 

sigxy:=  sqr(sigx)  +  sqr(sigy); 

radical:=  sqrt(sqr(sqr(sigx)  -  sqr(sigy))  +  4*sar(rho*sigx*sigy) ) ; 

sigS:=  sqrt((sigxy  -  radical)/2); 

sigL:=  sqrt((sigxy  +  radical)/2); 

cepmu:=  0.614*sigS  +  0.563*sigL; 

b:=  sqrt(sqr(mux)  +  sqr(muy)); 

v:=  b/cepmu; 

randc :=  cepmu*( 1.0039-0. 0528*v+0 . 4786*sqr ( v )-0 . 0793*atotheb(v , 3 ) ) 
end;  {RANDCEP} 


•**•*•******* 


{The  following  procedure  is  called  by  the  CBNORM 
procedure  which  calculates  the  CEPs  for  the 
correlated  bivariate  normal  (CBN)  approximation 
technique.  This  procedure  calculates  the  probability 
given  a  known  radius  (an  estimate  of  CEP).  It  is 
called  until  the  probability  is  within  0.0001  of 
0.5.  The  accuracy  of  this  method  can  be  improved 
by  increasing  the  number  of  intervals  used  to 
approximate  the  integral.  This  number  of  intervals 
is  a  user  input  that  is  asked  for  whenever  the  CBN 
method  is  used. 


procedure  CENORMPPOB ( var  sigx,  sigv,  mux,  muv,  rho,  bias:  real; 

ints:  integer); 

label  1; 

var  y,  hh,  i,  x,  c,  t,  a:  real; 
j  j  ,  j  ,  h :  integer; 

begin 

c :=  2*3. 1415927*sigy*sigx*sqrt( l-sqr(rho) ) ; 
c:=  1/c; 
i:=  bias; 
t  :=  0; 

for  j  :=  0  to  ints 
do  begin 

x:=  -i  +  ( j*(i/ ( ints/2) )) ; 
y:=  sqrt(abs(sqr(i)  -  sqr(x))); 
if  y  =  0  then  goto  1; 
for  h  :=  0  to  ints 
do  begin 

hh:=  -y  +  (h*(y/( ints/2) )) ; 

a:=  sqr( (x-mux)/sigx)-( 2*rho*( (x-mux)/sigx)*( (hh-muy )/sigy) )  + 
sq  r ( ( hh-muy ) / s igy ) ; 
a:=  a/( l-sqr( rho) ) ; 
if  a  >  140  then  a:=  0 
else  a:=  exp(-0.5*a); 
a : =  a*c; 

t:=  t  +  (a*( i/( ints/2) )*2*(y/ ints ) ) 
end ; 

1:  J3:-  3 

end ; 

cbnprob:=  t; 


********** 


{The  following  procedure  is  the  secant  method  applied 
to  the  correlated  bivariate  normal  approximation 
technique.  The  stopping  criterion  for  this  method 
is  when  the  probability  given  by  the  CBN0RMPR0B 
procedure  is  within  0.0001  of  0.5} 


procedure  CENORM  (var  sigx,  sigy,  mux,  muy,  rho,  cep: 

ints ) ; 


real ; 


var  cepipl,  cepl,  cepiml,  problpl,  probi,  probiml:  real; 


begin 

cepiml:=  cep; 

CBNORMPROB( s igx ,  sigy,  mux,  muy,  rho,  cepiml,  ints); 
probiml :=  cbnprob; 

if  cbnprob  <  0.5  then  cepi:=  1.05*cep 
else  cepi:=  0.95*cep; 

CBN0RMPR0B( s igx ,  sigy,  mux,  muy,  rho,  cepi,  ints); 
probi:=  cbnprob; 
probipl:=  0.5; 

while  abs(cbnprob  -  0.5)  >  0.0001  do 
begin 

cepipl :=  cepi+( ( cep iml-cep i )*( probipl -probi )/ (probiml -probi  ) ) ; 

CBNORMPROB(sigx,  sigy,  mux,  muy,  rho,  cepipl,  ints); 

cepiml:=  cepi; 

probiml :=  probi; 

cepi:=  cepipl; 

probi :=  cbnprob; 

end ; 

cbn:=  cepipl 
end ; 


********** 


{The  following  procedure  is  called  by  EXACT  which 

calculates  CEPs  for  the  "exact”  method.  EXACTPROB 
calculates  probabilities  given  an  estimate  of  the 
CEP.  It  is  called  until  the  probability  is  within 
0.0001  of  0.5.} 


procedure  EXACTPROB  (var  sigx,  sigy,  mux,  muy,  cep:  real); 

type  logarray  =  array [ 1 .. 100 ]  of  real; 

var  kcheck,  jcheck,  Lsum,  kf,  kmLf,  L2f,  hold,  La,  Lterm:  real; 
k2f,  kplf,  xka,  xk,  isum,  jf,  jmif,  i2f,  ia,  iterm:  real; 
j2f,  jplf,  yja,  yj ,  kpjplf,  jterm,  jsum,  ksum,  d:  real; 
r,  ja>  yjb,  xkb,  Lb,  ib,  Lc,  ic:  real; 
log:  logarray; 

i,  j,  k,  L,  il,  12,  jl,  j  2 ,  j  3 ,  j4,  kl,  k2,  k3,  LI,  L2 :  integer 


E-6 


.Va'V.V  •.  ’.’  V  V  '  '  O  -  • . 


begin  { EXACTPROB} 
r:  =  cep; 
kcheck:=  1; 
k:  =  0; 
ksum:®  0; 

while  kcheck  >  0.0001  do 
begin 
kf : =  0; 

if  k  >  0  then  for  kl:=  1  to  k  do  kf:=  kf  +  1 og [ k 1  ] ; 
k2f :=  0; 

if  k  >  0  then  for  k2:=  1  to  (2*k)  do  k2f:=  k2f  +  log[k2]; 
kplf:=  0; 

for  k3:=  1  to  (k+1)  do  kplf:=  kplf  +  log f k3 ] ; 
jcheck:=  1; 

j:=  0; 

jsum:=  0; 

while  jcheck  >  0.0001  do 
begin 

writeln(k,  '  ' ,  kcheck,  j,  '  jcheck); 

jf:=  0; 

if  j  >  0  then  for  jl:  =  1  to  j  do  jf:=  jf  +  log [.11]  ; 
jplf:=  0; 

for  j2:=  1  to  (j+1)  do  jplf:=  jplf  +  log[j2]; 
j  2  f :=  0; 

if  j  >  0  then  for  j3:=  1  to  (2*j)  do  j2f:*  j2f  +  loc[j3]; 
kpjplf:=  0; 

for  j4:  =  1  to  (k+j+1)  do  kpjplf:=  kpjplf  +  log[j4]; 

Lsum:=  0; 
for  L:=  0  to  k  do 
begin 

kmLf:*  0; 

if  (k-L)  >  0  then  for  Ll:=  1  to  (k-L)  do  kmLf:=  kmLf+log [LI ] ; 
L2f  :=  0; 

if  L  >  0  then  for  L2:=  1  to  (2*L)  do  L2f:=  L2f  +  log [ L2 ] ; 
hold:=  (2*sqr(mux) )/(sqr(sigx) ) ; 

La:»  atotheb(hold,  (L/2)); 

Lc : =  La; 

if  odd(L)  then  La:=  -La; 
if  (kf-kmLf-L2f)  <  -70  then  Lb:=  0 
else  Lb:=  exp(kf  -  kmLf  -  L2f); 
if  (abs(La)  <  1.0e-16)  or  (abs(Lb)  <  1.0e-16)  then  Lterir:=  0 
else  Lterm:=  La*Lb*Lc; 

Lsum:=  Lsum  +  Lterm 


( EXACT? ROE  continued} 

hold:=  sar(r )/ (8*( sq r( sigx) ) ) ; 
xka:=  3totheb(hold,  k); 
if  odd(k)  then  xka:=  -xka; 
if  (k2f-kplf-(2*kf ))  <  -70  then  xkb:=  0 
else  xkb:=  exp(k2f  -  kplf  -  (2*kf)); 
xk:=  xka*xkb*Lsum ; 
isum:=  0; 
for  i:=  0  to  j  do 
begin 

jmif:=  0; 

if  (j-i)  >  0  then  for  il:=  1  to  (j-i)  do  imif:=  jmif+log [ il ] ; 
i2f  :=  0; 

if  i  >  0  then  for  i2:=  1  to  (2*i)  do  i2f:=  i2f  +  log [ i2 ] ; 
hold :=  (2*(sqr(muy)))/(sar(sigy)); 
ia:=  atotheb(hold ,  (i/2)); 
ic : =  ia; 

if  odd(i)  then  ia:=  -ia; 
if  ( jf-jmif-i2f)  <  -70  then  ib:=  0 
else  ib:=  exp(jf  -  jmif  -  i 2 f ) ; 
if  (abs(ia)  <  1.0e-16)  or  (abs(ib)  <  1.0e-16)  then  iterm:=  0 
else  iterm:=  ia*ib*ic; 
isum:=  isum  +  iterm 
end ; 

hold : =  (sqr(r))/(8*(sqr(sigy))); 
yja:=  atotheb(hold ,  j); 
if  odd(i)  then  yja:  =  -yja; 
if  ( j  2f-j  plf-( 2*j f ) )  <  -70  then  yjb:*  0 
else  yjb : =  exp(j2f-  jplf  -  (2*jf)); 
yj:  =  yja*yjb*isurn ; 

if  (kplf+jplf-kpjplf )  <  -70  then  ja:=  0 
else  ja:  =  exp(kplf  +  jplf  -  kpjplf); 
jterm:=  xk*yj*ja; 
jsum:=  jsum  +  jterm; 

1check:=  abs(jterm); 

j:=  j  +  1 

end;  {while  jcheck  ...  } 
ksum:=  ksum  +  jsum; 
kcheck:=  abs(jsum); 
k :=  k  +  1 

end  {while  kcheck  ...  }  ; 

d :=  (exp(-(sqr(mux)/(2*sqr(s igx) ) )-( sqr(muy )/ ( 2*sqr(sigy ) ) ) ) )/ (sigx*sigy ) ; 
prob:=  (sqr(r)/2)*d*ksum 
end;  { EXACTPROB} 


********** 

{The  following  procedure  is  the  secant  method  applied 
to  the  "exact1’  approximation  technique.  It  calls 
the  procedure  EXACTPROB  until  the  probability  is 
within  0.0001  of  0.5.} 
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procedure  EXACT  (var  sigx,  sigy,  mux,  muv,  cep:  real); 


var  cepipl,  cepi,  cepiml,  probipl,  proM ,  probiml:  real; 


begin 

for  ii:=  1  to  100  do  log [ i i ] : =  ln(ii); 
cepiml:=  cep; 

EXACTPROB(sigx,  sigy,  mux,  muy,  cepiml); 
probiml : =  prob ; 

if  prob  <  0.5  then  cepi:=  1.05*cep 
else  cepi:=  0.95*cep; 

EXACTPROB(sigx,  sigy,  mux,  muy,  cepi); 
probi:=  prob; 
probipl:=  0.5; 

while  abs(prob  -  0.5)  >  0.0001  do 
begin 

cepipl :=  cepi+( (cepiml -cepi )*(probipl-probi )/ (probiml-probi ) ) ; 
EXACTPROB(sigx,  sigy,  mux,  muy,  cepipl); 
cepiml :=  cepi ; 
probiml :=  probi ; 
cepi:=  cepipl; 
probi : =  prob; 
end ; 

cep:=  cepipl 
end ; 
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{The  following  procedure  is  the  decision  maker  for  which 
of  the  algebraic  approximations  is  the  "best"  for  the 
given  regime  of  the  elliptic! ty/bias  parameter  space.} 


procedure  BESTALGEERAIC( var  ints:  integer; 

sigx,  sigy,  mux,  muy,  rho,  bias,  ellip:  real); 

begin 

if  ellip  <  0.05  then  writeln('SORRY,  YOU  MUST  CHOOSE  A  METHOD.'); 
if  (ellip  <  0.15)  and  (ellip  >=  0.05) 
then  begin 

if  bias  <=  200 
then  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 

write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 
wri teln( randc : 10 : 3 ) 
end ; 

if  bias  >  200  then  writeln('SORRY,  YOU  MUST  CHOOSE  A  METHOD.') 
end ; 

if  (ellip  <  0.25)  and  (ellip  >=  0.15) 
then  beg^n 

if  (bias  <  200)  or  ((bias  >=  350)  and  (bias  <  500)) 
then  begin 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 

write('THE  GRUBBS-PATNA IK/ CHI- SQUARE  IS  BEST  HERE:  '); 
writeln(gpcs : 10:3) 
end ; 

if  ((bias  >=  200)  and  (bias  <  350))  or 
((bias  >=  500)  and  (bias  <=  650)) 
then  begin 

RANDCEP( s igx ,  sigy,  mux,  muy); 
write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 
writeln(randc: 10:3) 
end; 

if  bias  >  650  then  writeln('SORRY,  YOU  MUST  CHOOSE  A  METHOD.') 
end ; 

if  (ellip  <  0.5)  and  (ellip  >=  0.25)  and  (bias  <  200) 
then  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 

wri teln( 'THE  MODIFIED  RAND-234  IS  EEST  HERE:  ',  randc:10:3) 
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{ BESTALGEERAIC  continued} 


if  (ellip  <  0.4)  and  (ellip  >=  0.25) 
then  begin 

if  ((bias  >=  200)  and  (bias  <  350))  or 

((bias  >=  500)  and  (bias  <  650))  or 

((bias  >=  900)  and  (bias  <  1000))  or 
((bias  >=  1200)  and  (bias  <  1300)) 
then  begin 

GPCSCEP( s igx,  sigy,  mux,  muy,  rho); 
write( 'THE  GRUEBS-PATNAIK/CHI-SQUARE  IS  BEST  HERE 
writeln(gpcs : 10:3) 
end ; 

if  ((bias  >=  350)  and  (bias  <  500))  or 

((bias  >-  650)  and  (bias  <  900))  or 

((bias  >=  1000)  and  (bias  <  1200))  or 
((bias  >=  1300)  and  (bias  <=  1550)) 
then  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 

write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 
writeln( randc : 10 : 3) 
end 

end ; 

if  (ellip  >=  0.4)  and  (ellip  <  0.5) 
then  begin 

if  ((bias  >=  200)  and  (bias  <  500))  or 

((bias  >=  600)  and  (bias  <  900))  or 

((bias  >=  1450)  and  (bias  <=  1550)) 
then  begin 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 
write( 'THE  GRUBBS-PATNAIK/CHI-SOUARE  IS  BEST  HERE 
wri teln (gpcs : 10 : 3 ) 
end ; 

if  ((bias  >=  500)  and  (bias  <  600))  or 

((bias  >=  900)  and  (bias  <  1450)) 

then  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 

write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 
wri teln ( randc : 10 : 3 ) 
end 

end ; 

if  (ellip  >=  0.5)  and  (ellip  <  0.7)  and 

((bias  <  600)  or  ((bias  >=  1300)  and  (bias  <  1450))) 
then  begin 

GPCSCEP( s igx ,  sigy,  mux,  muy,  rho); 

write( 'THE  GRUBBS-PATNAIK/CHI-SOUARE  IS  EEST  HERE:  ' 
writeln(gpcs : 10:3) 


{ BESTALGEERAIC  continued} 


if  (ellip  >=  0.5)  and  (ellip  <  0.6) 
then  begin 

if  ((bias  >=  600)  and  (bias  <  750))  or 
((bias  >=  1200)  and  (bias  <  1300)) 
then  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 

write('THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 
writeln( randc : 10 : 3) 
end ; 

if  ((bias  >=  750)  and  (bias  <  1200))  or 
((bias  >=  1450)  and  (bias  <  1550)) 
then  begin 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 

write( 'THE  GRUBBS-PATNAIK/CHI-SQUARE  IS  BEST  HERE:  ' 
writeln(gpcs : 10 : 3) 
end 

end ; 

if  (ellip  >=  0.6)  and  (ellip  <  0.7) 
then  begin 

If  (bias  >=  750)  and  (bias  <  900) 
then  begin 

RANDCEP( sigx,  sigy,  mux,  muy); 
write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 
writeln( randc: 10:3) 
end ; 

if  (bias  >=  950)  and  (bias  <  1300) 
then  begin 

GPWHCEP(sigx,  sigy,  mux,  muy); 

write('THE  GRUBBS-PATNAIK/WILSON-HILFERTY  IS  BEST  HERE 
wri teln(gpwh : 10 : 3) 
end 

end; 

If  (ellip  >=  0.6)  and  (ellip  <  0.65)  and 
(((bias  >=  900)  and  (bias  <  950))  or 
((bias  >=  1450)  and  (bias  <  1550))) 
then  begin 

GPWHCEP( s igx ,  sigy,  mux,  muy); 

wri te( 'THE  GRUBBS-PATNAIK/WILSON-HILFERTY  IS  BEST  HERE: 
writeln(gpwh: 10:3) 
end ; 

if  (ellip  >=  0.65)  and  (ellip  <  0.7)  and 
(((bias  >=  900)  and  (bias  <  950))  or 
((bias  >=  1450)  and  (bias  <  1550))) 
then  begin 

RANDCEP( s Igx ,  sigy,  mux,  muy); 

write('THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 

writeln(randc:10:3) 


{ BESTALGEBRAIC  continued} 


if  (ellip  >=  0.7)  and  (ellip  <  0.8) 
then  begin 

if  (bias  <  100)  or 

((bias  >=  950)  and  (bias  <  1050))  or 
((bias  >=  1450)  and  (bias  <  1550)) 
then  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 
write( 'THE  MODIFIED  RAND-234  IS  EEST  HERE:  '); 
write In ( randc : 10 : 3 ) 
end ; 

if  ((bias  >=  100)  and  (bias  <  950))  or 
((bias  >=  1300)  and  (bias  <  1450)) 
then  begin 

GPCSCEP( s igx ,  sigy,  mux,  muy,  rho); 

wri te( 'THE  GRUBBS-PATNA IK/ CHI-SQUARE  IS  BEST  HERE:  ') 
writeln(gpcs : 10 : 3) 
end ; 

if  (bias  >=  1050)  and  (bias  <  1300) 
then  begin 

GPWHCEP(sigx,  sigy,  mux,  muy); 

write('THE  GRUBBS-PATNAIK/WILSON-HILFERTY  IS  BEST  HERE 
writeln(gpwh: 10:3) 
end 

end ; 

if  (ellip  >=  0.8)  and  (ellip  <  0.9) 
then  begin 

if  ((bias  >=  100)  and  (bias  <  600))  or 
((bias  >=  750)  and  (bias  <  1050))  or 
((bias  >=  1300)  and  (bias  <  1450)) 
then  begin 

GPCSCEP( s igx ,  sigy,  mux,  muy,  rho); 

write('THE  GRUBBS-PATNA IK /CHI-SQUARE  IS  BEST  HERE:  ') 
writeln(gpcs : 10:3) 
end ; 

if  (bias  >=  1200)  and  (bias  <  1300) 
then  begin 

GPWHCEP(sigx,  sigy,  mux,  muy); 

write('THE  GRUBBS-PATNAIK/WILSON-HILFERTY  IS  BEST  HERE 
writeln(gpwh : 10:3) 
end ; 

if  (bias  >=  1450)  and  (bias  <  1550) 
then  begin 

RANDCEP( s igx ,  sigy,  mux,  nuy); 

write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 
wri teln ( randc : 10 : 3 ) 


{ EEST ALGEBRAIC  continued} 

if  (ellip  >=  0.8)  and  (ellip  <  0.85) 
then  begin 

if  (bias  <  100) 
then  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 

write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 

writeln(randc: 10:3) 
end ; 

if  ((bias  >=  600)  and  (bias  <  950))  or 
((bias  >=  1050)  and  (bias  <  1200)) 
then  begin 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 

write( 'THE  GRUBBS-PATNAIK/ CHI- SQUARE  IS  BEST  HERE:  '); 

writeln(gpcs : 10:3) 
end 

end ; 

if  (ellip  >=  0.85)  and  (ellip  <  0.9) 
then  begin 

if  (bias  <  100) 
then  begin 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 

write('THE  GRUBBS-PATNAIK/ CHI-SQUARE  IS  BEST  HERE:  '); 
writeln(gpcs : 10:3) 
end ; 

if  ((bias  >=  600)  and  (bias  <  750))  or 
((bias  >-  1050)  and  (bias  <  1200)) 
then  begin 

RANDCEP( s igx ,  sigy,  mux,  muy); 

write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 

write Ini randc: 10:3) 
end 

end ; 

if  (ellip  >=  0.9)  and  (ellip  <  0.95) 
then  begin 

if  (bias  <  600)  or 

((bias  >=  900)  and  (bias  <  1050))  or 
((bias  >=  1250)  and  (bias  <  1450)) 
then  begin 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 

write( 'THE  GRUBBS-PATNAIK /CHI-SOU ARE  IS  BEST  HERE:  '); 

writeln(gpcs : 10 : 3) 
end  ; 

if  ((bias  >=  600)  and  (bias  <  900))  or 
((bias  >=  1050)  and  (bias  <  1250))  or 
((bias  >=  1450)  and  (bias  <  1550)) 
then  begin 

RANDCEP( s igx ,  sigy,  mux,  muy); 

wri te( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 

writeln(randc: 10:3) 
end 

end ; 
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{ BESTALGEBRAIC  continued} 


if  (ellip  >=  0.95)  and  (ellip  <=  1.0) 
then  begin 

if  (bias  <  100)  or 

((bias  >=  200)  and  (bias  <  500))  or 

((bias  >=  750)  and  (bias  <  1200))  or 

((bias  >=  1300)  and  (bias  <  1450)) 
then  begin 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 

write( 'THE  GRUBBS-PATNAIK / CHI-SQUARE  IS  BEST  HERE:  ') 

wri teln (gpcs : 10 : 3) 
end ; 

if  ((bias  >=  100)  and  (bias  <  200))  or 

((bias  >=  500)  and  (bias  <  750))  or 

((bias  >=  1450)  and  (bias  <  1550)) 
then  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 

write( 'THE  MODIFIED  RAND-234  IS  BEST  HERE:  '); 

wri teln (randc : 10 : 3) 
end ; 

if  (bias  >=  1200)  and  (bias  <  1300) 
then  begin 

GPWHCEP(sigx,  sigy,  mux,  muy); 

write('THE  GRUBBS-PATNAIK/WILSON-HILFERTY  IS  BEST  HERE 

writeln(gpwh: 10:3) 
end 

end 

end;  {BESTALGEBRAIC} 


********** 


{The  following  procedure  is  used  to  execute  the  desired 
operation  that  the  user  decided  to  accomplish} 


procedure  PICKMETHOD(var  pick:  integer; 

sigx,  sigy,  mux,  muy,  rhc,  bias,  ellip:  real); 

beqin 

case  pick  of 
1:  begin 

writeln('HOW  MANY  INTERVALS  DO  YOU  VJANT  THF.  FOR  THE  CBN?'); 
writeln( 'THE  MORE  INTERVALS,  THE  MORE  TIME,  THE  MORE  ACCURATE'); 
readln( ints ) ; 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 

GPWHCEP(sigx,  sigy,  mux,  muy); 

RANDCEP(sigx,  sigy,  mux,  muy); 
cep:®  gpcs; 

CBNORM(sigx,  sigy,  mux,  muy,  rho,  cep,  ints); 

EXACT(sigx,  sigy,  mux,  muy,  cep); 

writeln('CEPs  FOR  ALL  METHODS  ARE  AS  FOLLOWS:'); 

writeln; 

wri te In ('GRUEES-PATNAIK/ CHI-SQUARE :  ',  gpcs: 10:3); 

wri teln( 'GRUBES-PATNAIK/WILSON-HILFERTY :  ',  gpwh:10:3): 

write In ('MODIFIED  RAND-234;  ',  randc:10:3); 

writeln('CEN:  ',  cbn:10:3); 

wri teln( '"EXACT" :  ',  cep:10:3) 

end ; 

2:  begin 

GPCSCEP(sigx,  sigy,  mux,  muy,  rho); 

write( 'THE  CEP  GIVEN  BY  THE  GRUEBS-PATNAIK/ CHI-SOUARE  IS:  '); 
writeln (gpcs: 10:3) 
end ; 

3:  begin 

GPWHCEP(sigx,  sigy,  mux,  muy); 

wri te( 'THE  CEP  GIVEN  BY  THE  GRUBBS-PATNAIK/WILSON-HILFERTY  IS:  '); 
writeln(gpwh: 10:3) 
end; 

4:  begin 

RANDCEP(sigx,  sigy,  mux,  muy); 

writeln('THE  CEP  GIVEN  BY  THE  MODIFIED  RAND-234  IS:  ',  randc:10:3); 
end ; 

5:  begin 

writeln('HOW  MANY  INTERVALS  DO  YOU  WANT  THE  FOR  THE  CBN?'); 
writeln( 'THE  MORE  INTERVALS,  THE  MORE  TIME,  THE  MORE  ACCURATE'); 
readln( ints ) ; 

CBNORM(sigx,  sigy,  mux,  muy,  rho,  cep,  ints); 
writeln( 'THE  CEP  GIVEN  BY  THE  CBN  IS:  ',  cbn:10:3) 
end ; 

6 :  beg i n 

EXACT(sigx,  sigy,  mux,  muy,  cep); 

writeln( 'THE  CEP  GIVEN  BY  THE  "EXACT"  IS:  ',  cep: 10:3) 
end ; 

7:  BESTALGEBRAIC( ints ,  sigx,  sigy,  mux,  muy,  rho,  bias,  ellip) 
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{The  following  is  the  Main  program.  It  accepts  the  input 
of  the  user  for  the  variables  and  then  which  method 
the  user  would  like  to  calculate} 

begin  {MAIN  PROGRAM) 

writeln( 'PLEASE  INPUT  THE  SAMPLE  STANDARD  DEVIATIONS  AND  MEANS'); 
write( 'CROSSRANGE  MEAN  {X  BAR}:  '); 
readln(mux) ; 

write( 'CROSSRANGE  STANDARD  DEVIATION  {SIGMA  X} :  '); 

readln(sigx) ; 

wri te( 'DOWNRANGE  MEAN  {Y  BAR}:  '); 
readln(muy) ; 

write ('DOWNRANGE  STANDARD  DEVIATION  {SIGMA  Y} :  '); 

readln(sigy) ; 

write('INPUT  THE  CORRELATION  COEFFICIENT,  IF  ANY  {RHO}:  '); 
readln(rho) ; 

bias:=  sqrt(sqr(mux)  +  sqr(muy)); 
ellip:=  sigx/sigy; 

writeln('HOW  WOULD  YOU  LIKE  FOR  THE  PROGRAM  TO  RUN?'); 
wri teln( '  {ENTER  APPROPRIATE  NUMBER}'); 

writeln; 

writeln('l  =  CALCULATE  ALL  METHODS'); 

writeln('2  =  CALCULATE  GRUBBS-PATNAIK/ CHI-SQUARE  METHOD  ONLY'); 

wri teln('3  =  CALCULATE  GRUBBS-PATNAIK/ WILSON-HILFERTY  METHOD  ONLY'); 

wri teln( '4  *  CALCULATE  MODIFIED  RAND-234  METHOD  ONLY'); 

writeln('5  -  CALCULATE  CBN  METHOD  ONLY'); 

writeln('6  =  CALCULATE  "EXACT”  METHOD  ONLY'); 

writeln('7  =  CALCULATE  THE  "BEST"  OF  THE  ALGEBRAIC  METHODS'); 

readln(pick) ; 

writeln; 

writeln( "*************■***************★*****■**★**********★**************' )  j 
writeln ; 

PICKMETHOD(pick,  sigx,  sigy,  mux,  muy,  rho,  bias,  ellip); 
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